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ABSTRACT 

It  is  known  that  the  static  and  dynamic  behavior  of 
fluids  in  porous  media  depends  to  a  large  measure  on  porous- 
media  geometry.   In  the  past,. the  ability  to  characterize 
this  geometry  has  been  restricted  to  such  average  properties 
as  porosity  and  permeability.   However,  in  recent  years 
attempts  have  been  made  to  achieve  a  more  precise  charac- . 
terization  based  upon  the  fact  that  porous  media  is 
statistically  composed. 

In  this  thesis  techniques  from  statistical  'communi-' 
cation  theory  are  adapted  as  possible  methods  for  accomplish' 
ing  this  classification  process.   In  simplest  form,  a 
dichotamous  function  is  defined  by  passing  a.  line  through  a 
porous  medium,  the  function  having  one  value  when  the'  line 
is  in  solid  matrix,  and  another  value  when  the  line  passes 
through  pore  space.   The  function  is  then  analyzed  using 
(1)  Classical  Fourier  series  harmonic  analysis,  and  (2) 
determination  of  the  autocovariance  estimate  and  power 
spectrum.   Comparisons  are  made  between  seyeral  functions 
created  from  the  same  medium,  and  with  functions  created 
from  other  media. 

The  results  indicate  that  the  autocovariance  estimate  * 
and  the  power  spectrum,  as  characterizing  functions,  can 
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discriminate  between  different  media.   This  success  suggests 
many  more  paths  of  investigation,  possibly  leading  to  the 
complete  characterization  of  porous  media  through  statistical 
analysis . 
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CHAPTER  I 

INTRODUCTION 

In  August,  1965»  two  American  astronauts  journeyed 
over  two  hundred  miles  into  space  in  preparation  for  man's 
future  journeys  to  the  moon  and  thousands  of  miles  beyond. 
As  of  August,  1965)  man  had  journeyed  only  a  few  hundred 
feet  beneath  the  surface  of  the  earth  upon  which  he  lives. 
In  this  very  interesting  paradox,  man  is  now  capable  of 
gathering  information  many  miles  into  space,  but  has  very 
little  first-hand  information  about  the  sub-surface  of  the 
very  earth  from  which  he  launches  his  space  probes.   If  this 
seems  incongrous,  consider  that  while  the  benefits  of  space 
exploration  are  speculative  at  best,  man  derives  the  vast 
majority  of  his  resources,  his  energy  sources,  and  even  his 
food,  from  this  earth.   However,  there  are  many  explanations 
for  man's  actions.   Perhaps  the  foremost  is  the  extreme' 
difficulty  encountered  in  any  attempt  to  explore  extensively 
beneath  the  earth's  surface.   Costs  would  certainly  be  pro- 
hibitive, even  if  the  engineering  difficulties  could  be 
overcome.   Another  explanation  might  be  that  man  has  been 
able  to  theorize  sufficiently  about  the  earth's  sub-surface 
areas  from  the  limited  information  available  so  that  he  has 
not  been  willing  or  seen  the  necessity  of  risking  great 
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expense,  for  what  might  be  very  little  gain.   A  less  prac- 
tical explanation,  but  in  all  probability  a  very  real  reason 
for  the  lack  of  extensive  sub-surface  exploration,  would  be 
that  there  is  less  glamour  to  most  people  about  traveling  to 
the  center  of  the  earth  than  in  a  journey  to  Venus  or  Mars 
and  thus,  less  interest  in  making  the  attempt. 

The  above  discussion  is  not  meant  as  an  implication  that 
man  is  not  interested  in  the  earth's  sub-surface.   To  the 
contrary,  geologists,  petroleum  engineers,  and  others  who 
derive  their  livelihood  from  the  soil  and  beneath  the  soil, 
have  performed  extensive  research  in  this  area.   From  the 
exploration  of  caves  and  mines,  from  core  samples  drawn  from 
wells,  and  in  many  other  ways,  man  has  been  able  to  learn 
much  about  the  earth's  constitution  and  composition.   As  a 
result  of  this  research,  great  improvements  have  been  realized 
in,  the  methods  used  to  extract  the  natural  resources  from 
the  earth,  not  to  mention  the  improved  conservation  procedures 
which  will  allow  man  to  benefit  from  the  earth's  resources 
for  many  years  to  come. 

As  a  part  of  this  research,  man  has  spent  many  years 
developing  intelligence  concerning  the  physical  properties 
of  porous  media  found  in  the  strata  of  various  kinds  that 
lie  beneath  the  earth's  surface  and  the  behavior  of  fluids 
within  these  porous  media.   Research  efforts  have  been  con- 
centrated on  four  properties  which  are  considered  fundamental 
and  the  basis  for  discussions  regarding  all  other  properties. 
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These  four  properties  are  porosity,  a  measure  of  the  void 
space  in  the  media;  permeability,  a  measure  of  the  media's 
ability  to  permit  flow;  the  fluid  saturation,  and  the 
electrical  conductivity  of  the  media.   We  will  discuss  these 
properties  at  greater  length  in  Chapter  II. 

The  basis  for  this  report  is  derived  from  a  paper  by 
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A.  E.  Scheidegger   on  "Statistical  Geometry  of  Porous  Media," 

published  in  1961.   As  pointed  out  by  Scheidegger,  although 
the  usual  geometrical  quantities,  such  as  porosity  and 
permeability,  give  generally  good  descriptions,  there  is  a 
need  for  geometrical  description  of  greater  detail.   This 
need  arises  because  of  the  extreme  complexity  of  the  subject 
matter.   For  example,  in  some  cases,  media  with  a  given 
porosity  will  hot  allow  as  much  flow  as  other  media  with  less 
void  space  because  of  the  lack  of  continuity  in  the  pore 
spaces  of  the  first  sample.   Media  with  identical  permea- 
bilities to  one  fluid  may  have  different  permeabilities  to 
another  fluid.   These  are  well-known  phenomena  in  the  area 
of  petroleum  engineering.   As'  a  result,  it  is  sometimes  very 
difficult  to  classify  media  of  different  types  without 
relying  on  extensive  laboratory  tests,  which  are  expensive 
and  time  consuming.   It  is  the  aim  of  this  report,  to  take 
advantage  of  .the  inherent  statistical  nature  of  porous  media 
to  characterize  the  media  by  a  set  of  numbers  or  a  mathe- 
matical function.   To  clarify  the  above  statement,  it  is 
important  to  point  out  that  porous  media  are  consolidated  in 


4 
a  statistical  nature.   An  individual  particle  has  properties 
determined  by  its  size,  shape,  density,  electrical  conduc- 
tivity, and  various  other  characteristics.   The  total  media 
will  have  properties  determined  by  the  properties  of  the 
individual  particles,  and  by  the  means  in  which  the  individual 
particles' have  united  to  form  the  whole.   This  consolidation 
process  will  in  turn  depend  upon  the  arrangement  of  the 
elementary  particles,  the  relative  number  of  different  types 
of  particles,  and  upon  the  relative  number  of  large  and- 
small  particles.   This  dependency  upon  distribution  functions 
is  obviously  statistical  in  nature. 

While  much  time  has  been  devoted  to  research  on  the 
various  physical  properties  mentioned  above,  very  little  has 
been  accomplished  in  the  latter  area  of  investigating  the 
statistical  nature  of  porous  media.   Scheideggar  and  Flood 
have  provided  impetus  in  this  area,  and  defined  some  of  the 
goals.   It  would  be  particularly  desirable,  for  example,  to 
be  able  to  look  at  two  pieces  of  porous  media  through  a 
microscope  and  tell  if  they  are  identical  in  a  statistical 
sense.   It  would  be  even  more  desirable  to  be  able  to  make 
intelligent  predictions  about  the  other  fundamental  properties 
from  the  statistical  description. 

The  purpose  of  this  report  is  to  introduce  and  attempt 
techniques  which  utilize  statistical  communication  theory  to 
characterize  porous  media.   In  simplest  form,  the  method  will 
consist  of  defining  a  function  created  by  passing  a  line  over 
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or  through  a  porous  medium,  and  assigning  a  value  of  +l,when 
the  line  passes  through  pore  space,  and  -1  when  it  passes 
through  solid  matter.   This  function  will  then  be  analyzed 
in  a  number  of  ways.   One  attempt  will  be  to  fit  the  function 
with  a  single  Fourier  series  and  determine  whether  the 
Fourier  coefficients  effectively  characterize  the  medium. 
Brief  mention  will  be  made  of  the  double  Fourier  series 
whereby  a  surface  area,  rather  than  a  single  line,  might  be 
analyzed. 

Finally,  the  function  will  be  looked  at  in  terms  of  cer- 
tain statistical  properties,  derived  from  communication 
theory.   The  idea  of  random  functions  will  be  introduced  and 
characterizing  functions  such  as  the'  autocovariance  function 
and  the  power  spectrum  will  be  determined.   Comparisons  will 
be  made  within  the  same  medium,  and  against  different  media. 
From  the  results,  it  is  hoped  to  ascertain  whether  a  given 
porous  medium  can  be  described  and  effectively  classified  by 
its  power  spectrum  and  autocovariance  function,  to  the' 

exclusion  of  all  others. 
i- 

If  any  of  these  three  methods  can  be  used  to  effectively 
classify  porous  media,  the  obvious  step  to  be  taken  later 
will  be  the  correlation  of  these  descriptive  quantities  with 
the  fundamental  properties,  such  as  porosity  and  permea- 
bility, which  are  now  used  to  characterize  porous  media. 


CHAPTER  II 

HISTORY 

Until  quite  recently,  man  has  investigated  and  classi- 
fied porous  media  largely  according  to  four  fundamental 
properties  from  which  it  was  felt  all  other  important  prop- 
erties could  be  derived.   These  four  basic  properties  are 
permeability,  porosity,  fluid 'saturation ,  and  electrical 
conductivity.   Some  of  the  earliest  work  was  done  by  Henry 
Darcy,  a  Frenchman.   In  I856  Darcy ,  investigating  the  flow 
of  water  through  a  sand  filter,  equated  flow  rate  with  the 
cross  sectional  area  of  the  conduit,  the  length,  and  the 
pressure  head  created  by  having  the  ends  of  the  conduit  at 
different  levels.   He  found  that  a  proportionality  constant 
was  required  to  balance  the  equality  for  various  filters,  and 
that  the  constant  was  a  function  of  the  sand  packing.   This 
constant  became  known  as  the  permeability  of  the  media.   The 
unit,  of  permeability  is  appropriately  called  the  darcy,  and 
is  defined  in  the  petroleum   industry  as: 

A  porous  medium  has  a  permeability  of  one  darcy 
when  a  single  phase  fluid  of  one  centipoise  viscosity 
that  completely  fills  the  voids  of  the  medium  will 
flow  through  it  under  conditions  of  viscous  flow  at  a 
rate  of  one  cubic  centimeter  per  second  per  square 
centimeter  cross  sectional  area  under  a  pressure  or 
equivalent  hydraulic  gradient  of  one  atmosphere  per 
centimeter . 
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As  can  be  seen  from  the  definition,  Darcy's  law  may  be 
extended  to  include  other  fluids  as  well  as  water.   Ideally, 
the  permeability  is  a  function  only  of  the  porous  medium 
and  not  the  particular  fluid  flowing.   However,  investigators 
such  as  Klinkenberg  reported  permeability  variations  dis- 
covered when  using  gases  as  the  flowing  fluid.   The  permea- 
bility was  found  to  vary  with  the  gas  pressure,  and  the 
variation  was  attributed  to  a  phenomenon  known  as  gas  slip- 
page.  A  discussion  of  this  phenomenon  is  presented  in  most 
standard  texts,  e.g.,  Amyx,  e_t  al  ♦  ' 

Prior  to  discussing  the  work  of  Carman  and  Kozeny  in  the 
prediction  of  the  permeability  of  a  porous  media  from  its 
known  pore  size  distribution,  we  should  discuss  early  work  in 

the  area  of  calculating  porosity.   The  early  investigations  of 
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Slichter  and  Fraser  and  Graton   were  concerned  with  the 

prediction  of  porosity  of  packings  of  uniform  spheres.   Since' 

porous  media  are,  in  most  cases,  quite  obviously  not  composed 

of  uniform  spheres,  but  of  a  myriad  of  particle  sizes  and 

shapes  consolidated  with  every  conceivable  angularity,  the 

work  of  these  men  must  be  enlarged  upon  to  deal  correctly  with 

real  materials.   Whereas  their  work  was  largely  theoretical, 

and  others  continued  to  work  with  predictions  of  porosity  based 

upon  grain  size  distribution  utilizing  seive  analysis,  the 

more  recent  work  in  this  area  has  been  through  empirical 

determinations  of  porosity,  by  laboratory  measurement,  using 

such  devices  as  the  Kobe  porosimeter    (mercury  introduction), 
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or  the  Washburn-Bunting   porosimeter  (air  extracted  by 
vacuum  and  measured). 

Studies  of  pore  size  distribution  have  been  carried  out 
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by  several  investigators  including  Carman,    using  gas 

13 
adsorption  techniques;  Ishkin  and  Kazaner,    passing  two 

1^ 
phase  fluids  through  the  media;  Hitter  and  Erich,    using 

small  angle  x-ray  scattering  patterns,  and  by  many  others 
with  variety  of  techniques.   Grain  size  distribution  work 
has  been  accomplished  mostly  by  simple  seive  analysis  and 
there  is  apparently  a  need  for  more  work  in  this  area.   In 
fact,  Nasson  ^  indicates  that  after  extensively  researching 
this  field,  there  is  "a  definite  need  for  an  inclusive 
mathematical  function  containing  a  few  definite  parameters 
that  would  define  a  given  size  distribution  for  most  con- 
ditions approaching  ideal  grinding." 

One  of  the  most  interesting  works  done  in  the  deter- 
mination of  pore  size  distribution  is  attributed  to  Nuss  and 

-i  z 

Whiting,  ,  who  devised  a  method  of  impregnating  core  samples 
with  plastic  under  high  pressure.   The  plastic  was  then 
polymerized  at  a  high  temperature  and  the  rock  material 
removed  by  acid  leaching.   The  remaining  plastic  reproduces 
the  original  pore  space. 

Many  investigators  have  completed  correlations  between 

17 
the  porosity  and  permeability  of  porous  media.   Carman, 

Kozeny,   and  Wyllie   developed  predictive  equations  from 

these  correlations,  whereby  one  property  might  be  estimated 


accurately  from  the  other.   In  recent  years  in  England, 

20 
Collls-George  and  Childs   have  refined  the  Kozeny  Equation 

to  improve  results  and  allow  greater  generality.   Francher, 

21 

Lewis  and  Barnes    experimented  with  various  porous  matter 

to  determine  the  relationship  between  grain  size  and  fluid 
conductance.   In  their  work  they  accounted  for  both  laminar 

and  turbulent  flow.   A  brief  discussion  of  the  results  is 

22 
found  in  the  text  by  Amyx . 

Nearly  all  rock  samples  brought  to  the  surface  contain 
amounts  of  liquid  entrained  in  their  pore  space.   The  deter- 
mination of  the  amounts  and  types  of  liquid  is  a  very 
important  characteristic  description  of  the  porous  media. 

Comprehensive  descriptions  of  the  methods  used  to  accomplish 

2  3        2k 
this  determination  are  covered  by  Hill,  ^  Horner, 

25 
Schilthius,    and  others.   The  two  most  common  techniques 

for  determination  of  fluid  content  are  by  retorting,  and  by 
solvent  extraction  and  distillation.   Fluid  saturation  may 
also  be  determined  by  using  an  indirect  approach,  that  is, 
by  measuring  some  other  property  and  correlating  the  findings 
with  past  data  to  estimate  the  fluid  entrained  within  the 
sample.   Use  of  electric  logs  and  capillary  pressure  measure- 
ments are  two  such  indirect  means  of  arriving  at  the  fluid 
saturation. 

The  fourth  fundamental  property  of  porous  rock  is  the 
electrical  conductivity.   The  rock,  as  it  comes  from  the 
ground,  is  comprised  of  solid  fragments,  void  space  and  some 
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amount  of  entrained  liquid.   Generally,  the  solids  and  the 
oil  and  gas  components  filling  the  pore  spaces  are  non- 
conductors of  electricity.   The  water,  if  it  contains  dis- 
solved salts,  then,  is  the  only  conductor.   The  ability  or 
non-ability  of  a  rock  to  conduct  electricity  is  thus  an 

indication  of  the  interconnected  or  continuous  water  filled 

26 
pores.   Wyllie  and  Spangler   created  the  first  model  of 

porous  media  taking  advantage  of  this  phenomenon  to  charac- 
terize the  media  according  to  its  electrical  resistivity, 
and  hence,  derived  such  other  properties  as  tortuosity  and 
porosity.   Cornell  and  Katz  '   and  Wyllie  and  Gardner 
added  refinements  in  later  models. 

From  the  four  basic  properties  covered,  many  other 
properties  can  be  derived.   Researchers  in  many  fields  have 
developed  various  descriptive  properties  which  characterize 
media  more  explicitly  for  their  particular  use.   For  example, 
the  geologist  has  expanded  the  characterization  to  include 
mineral  content,  while  the  petroleum  engineer,  because  of  his 
interest  in  flow  properties,  has  developed  such  concepts  as 
capillary  pressure,  relative  permeability,  and  wettability. 
A  description  of  these  properties  is  beyond  the  scope  of  this 
report,  but  the  bibliography  references  the  works  in  this 

area  of  Purcell,2^  Leverett,^0  H.  Brown, ^1  Rose,^2  Slobod,-^ 

3/4.  35 

Fatt,   and  Rapaport  and  Leas .  ^  % 

Man  has,  then,  taken  advantage  of  many  physical  properties 

to  descriptively  classify  porous  media.   However,  in  1961, 
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A.  E.  Scheidegger  and  H.  D.  Fara,    writing  in  the  Journal  of 

Geophysical  Research,  point  out  that  the  usual  geometrical 

properties,  such  as  porosity  and  permeability,  do  not  describe 

the  media  sufficiently  for  present  requirements.   Quoting 

their  article: 

Geometrically,  a  porous  medium  is  defined  by 
giving  the  analytical  equation  of  the  surface  which 
bounds  its  pore  space.   For  any  practical  purpose, 
this  is  impossible  to  accomplish.   Average  geometrical 
quantities,  such  as  porosity  .  .  .  ,  and  specific 
surface  .  .  .  ,  are  conceptually  acceptable  (and  can 
be  measured),  but  this  is  not  true  for  the  concept  of 
"pore  size."   Since  the  pores  form  an  extremely  compli- 
cated system  .  .  .  ,  it  is  not  easy  to  determine  a 
pore  size  in  a  conceptually  satisfactory  fashion. 

Nevertheless,  there  is  a  great  need  for  esta- 
blishing a  better  geometrical  characterization  of  a 
porous  medium  than  is  afforded  by  the  mere  deter- 
mination of  porosity  and  specific  surface  area. 

Scheidegger  and  Fara  desired  a  means  of  characterizing  porous 

media  through  analysis  of  the  statistical  make-up  of  the 

media. 

These  investigators  introduced  four  possible  methods  of 
statistical  analysis:   (1)  by  use  of  a  correlation  function, 
(2)  a  spectral  analysis  in  terms  of  harmonic  functions,  (3)  a 
spectral  analysis  in  terms  of  other  orthogonal  functions, 
and  (4)  a  spectral  analysis  of  a  specially  constructed  charac- 
teristic function  of  the  media. 

In  each  of  the  four  techniques,  the  writers  developed 
statistical  data  for  interpretation  by  drawing  an  arbitrary 
line  through  the  porous  medium.   The  line  has  a  designated 
origin  and  data  points  are  selected  at  equally  spaced  points 
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along" the  line.   At  each  of  these  points,  a  value  of  +1  is 
assigned  if  the  line  is,  at  that  point,  passing  through  void 
or  pore  space,  and  a  value  of  -1  is  assigned  if  the  line  is 
passing  through  solid.   This  series  of  +1  and  -1  values, 
taken  from  the  origin,  define  an  orthogonal  function  f(s), 
which  is  characteristic  of  the  line  drawn  through  the  function 
and  of  the  media  itself.   If  the  medium  is  homogeneous  and 
isotropic,  and  if  the  sample  length  is  of  sufficient  size, 
theoretically,  the  functions  or  lines  are  independent  of 
direction,  origin,  or  position.   The  function  is  treated  as 
a  random  function,  as  defined  by  Moyal,--5'  and  can  be  statis- 
tically analyzed  in  several  ways. 

If  the  function  is  defined  as  above  for  void  and  solid 
spaces,  the  mean  of  the  function  may  be  derived: 

f  ■  n(P)  +  (-D(i-p)  '  (2.D 

where  P  7  porosity  or  pore  spaces 

.'.  T   =  2P"1  (2.2) 

In  many  cases,  Scheidegger  points  out  that  it  is  incon- 

venient  to  work  with  a  function  having  a  mean  other  than  zero,  * 

so  he  rewrites  (2)  as  a  new  function  with  a  mean  that  is 

equal  to  zero.   The  new  function  becomes: 

f'(s)  =  f(s)  -  (2P-1)  (2.3) 

In  the  first  of  the  four  techniques'  which  the  authors 
■  » 

suggest,  the  function  is  characterized  by  its  autocorrelation 

function: 
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R(A)  =  lira    S*l   f(s)  f(s+A)ds  ,-  M 

S >°°  -=§ T \*»i*) 

;_s  ds 

Chapter  III  of  this  paper  contains  a  discussion  of  the  auto- 
correlation function. 

The  second  technique  referred  to  in  the  Scheidegger 
paper  is  to  fit  the  function,  f(s),  with  a  Fourier  series, 
and  evaluate  the  Fourier  coefficients,  a,  and  b,  .   A  suitable 
evaluation  suggested  is  to  use  the  equation: 

ck2  -  ak2  ♦  b/  (2.5) 

where  c,   is  not  only  a  fair  representation  of  the  medium's 
geometry  (statistically),  but  has  the  added  advantage  of  being 
independent  of  point  of  origin   or  phase  angle. 

The  third  and  fourth  techniques  suggested  are  mere 
refinements  of  the  two  discussed  briefly  above.   Each  of  the 
four  techniques  was  offered  by  Scheidegger  and  Fara  in  an 
introductory  manner,  and  there  was  very  little  in  the  way  of 
tangible  conclusions.   However,  these  men  opened  an  entirely 
new  facet  to  the  subject  of  characterizing  porous  media. 
Unfortunately,  in  correspondence  between  the  author  of  the 
report  in  hand,  and.  Dr.  Scheidegger,  the  latter  indicated  that 
he  had  done  no  further  work  in  this  area.   Further  works  in 
the  refinement  of  these  and  similar  techniques  are  now  left 
to  this  author  and  others. 

Perhaps  brief  mention  should  be  included  herein  as  to 
the  history  of  statistical  communication  theory  and  time 
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series  analysis,  which  Scheidegger  has  introduced  into  the 
world  of  porous  matter.   Granger   credits  the  astronomers 
as  being  the  first  to  analyze  time  series,  although  he  men- 
tions discussions  of  meteorological  and  geophysical  series 
in  the  middle  of  the  nineteenth  century.   Even  in  these  early 
discussions,  it  was  noted  that  series  of  data  taken  over 
periods  of  time  often  follow  patterns  and  contain  certain 
trends,  sometimes  obvious  and  sometimes  not  so  obvious.   If 
the  obvious  trends  such  as  those  resulting  from  secular  or 
seasonal  variations  could  be  subtracted  from  the  series,  the 
remaining  series  would  contain  random  and  independent  data, 
and  perhaps  some  of  the  not  so  obvious  trends  or  patterned 
variations . 

Many  of  the  works  in  separating  these  "hidden  periodicities" 

39  40  4l 

were  instituted  by  La  Grange,   Buys-Ballot,   Whittaker, 

42  43 

Stokes,   and  Schuster.  J     In  all  of  these  works,  the  random  . 

independent  series  was  assumed  to  be  unimportant,  and  attempts' 

were  made  to  mathematically  generate  similar  functions  as 

those  creating  the  hidden  trends. 

Schuster  initiated  the  use  of  the  periodogram,  a  method 

of  pinpointing  lesser  periodic  variations  in  a  series.   This 

approach  was  refined  by  such  investigators  as  Whittaker  and 

44 
Robinson   who  showed  that  apparent  major  trends  may  only  be 

a  combination  of  multiple  lesser  trends  at  their  coincident 

peaks . 

4S        46         '47 
It  followed  that  Kolmogoroff,  J   Weiner,    and  Cramer  ' 

should  refine  these  attempts  at  separating  periodic  tendencies 
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by  using  Fourier  series  analysis  to  segregate  obvious  sine 
or  cosine  trends.   All  of  these  earlier  theoretical  works 
led  to  practical  applications,  particularly  in  the  field  of 
communications  by  Tukey   and  his  assistants  at  the  Bell 
Telephone  Laboratories.   Applications  in  other  fields  have 
followed,  until  today,  time  series  analysis   and  communication 
theory  are  major  tools  of  the  business  and  economic  world 
as  well  as  the  world  of  science  and  engineering. 


CHAPTER  ITI 
HARMONIC  ANALYSIS  OF  RANDOM  FUNCTIONS 

A  simple  means  of  analyzing  a  function  is  through  a 
resolution  into  its  harmonics  by  the  use  of  a  Fourier 
series.   Where  the  function  has  one  independent  variable 
and  one  dependent  variable,  f(x)  =  z,  the  one  dimensional 
or  single  Fourier  series  is  appropriate.   A  two  dimen- 
sional function  with  two  independent  variables  and  the  one 
dependent  variable,  f(x,y)  =  z,  requires  the  use  of  the 
double  Fourier  series.   A  discussion  of  these  mathematical 
tools  is  contained  in  Appendixes  A  and  B.   Where  random 
processes  or  random  functions  are  being  considered,  the 
classical  harmonic  analysis  is  not  directly  applicable. 
However,  it  has  been  shown  that  a  "type"  of  Fourier  analysis 
is  still  possible.   This  theory,  which  has  .grown  out  of 
statistical  theory  of  communication,  provides  the  engineer 
with  additional  tools  to  assist  in  analyzing  random 
phenomena. 

Two  very  important  properties  of  a  signal  or  function, 
which  for  our  purposes  serve  to  characterize  the  function, 
are  the  power  spectrum  and  the  autocorrelation  function. 
Prior  to  presenting  the  mathematical  background  for  cal- 
culating these  properties,  it  may  be  well  to  offer  a  lay- 
men's definition  of  each. 

16 
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A.   Power  Spectrum 

This  is  a  term  used  in  basic  communication  theory  which 
has  been  borrowed  from  the  field  of  electricity.   Various 
analogies  have  been  drawn  to  describe  the  difficult  concept 
of  spectral  analysis,  one  of  the  better  ones  being  that 
discussed  in  Granger  and  Hatana.ka.     The  author  will 
attempt  to  draw  an  analogy  of  similar  nature.   If  a  person 
stands  on  the  corner  of  42nd  Street  and  Broadway,  during  a 
very  busy  period,  his  ears  will  pick  up  many  and  varied 
sounds.   Each  of  these  sounds  will  appear  to  be  different, 
although  in  many  cases  the  "sound"  which  one  hears  is  in 
reality  a  combination  of  many  sounds  coming  from  a  range  of 
signal  creating  devices. 


rat-. 
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Figure, 3«1   The  Human  Ear  Performing  Simple 
Spectal  Analysis 
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The  different  sounds  are  a  result  Df  many  signal 
generators,  emitting  different  frequencies.   The  human 
brain  is  equipped  with  a  device,  not  unlike  a  filter  used 
in  electrical  circuits,  which  isolates  each  signal  and 
attempts  to  analyze  it.   Sometimes  it  is  successful  and 
reports,'  "That  is  an  automobile  horn,  and  I  must  get  out  of 
the  street,"  and  sometimes  the  cacophony  is  so  great  that 
the  individual  signals  can  not  be  separated  and  Interpreted. 
The  brain  is  performing  a  type  of  spectral  analysis  when  it 
selects,  from  a  number  of  frequencies  that  the  ear  is  hear- 
ing, individual  signals  and  classifies  them  according  to 
amplitude  and  frequency.   Of  course,  in  communication  theory, 
and  perhaps  even  in  our  analogy,  there  may  be  so  many  fre- 
quencies that  it  is  only  within  our  capability  to  identify 
and  classify  within  ranges  of  frequencies.   This  is  more 
often  the  case  where  a  receiver  is  obtaining  both  intended 
or  predetermined  signals  and  "noise"  which  is  added  by  some 
unknown  and  random  generator.   The  mathematical  process  of 
decomposing  a  composite  signal  into  its  component  frequencies 
and  determining  the  proportion  of  the  total  signal  at  each 
frequency  has  been  called  Power  Spectrum  Analysis,  or  it 
is  often  abbreviated  as  Spectral  Analysis. 

Spectral  analysis  has  been  used  more  and  more  in  recent  . 
years  to  analyze  economic  time  series  for  the  business 
world.   Brown"5  describes  the  use  of  this  theory  to  separate 
the  many  factors  which  influence  an  economic  series,  such  as 
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seasonal  and  secular  trends.   Here,  as  in  electrical  and 
communication  work,  we  discover  the  presence  of  "noise," 
or  factors  which  we  cannot  interpret  as  to  origin,  or 
predict  as  to  future  behavior. 
B.   Autocorrelation  Function 

The  autocorrelation  function  may  be  described  as  the. 
new  function  created  when  a  function  is  correlated  with 
itself.   However,  in  order  to  provide  a  more  meaningful  and 
definitive  picture,  we  must  look  at  the  definition  shown  in 
Ekeblad .    This  text  points  out  that  the  degree  of  auto- 
correlation in 'a  time  series  indicates  how  well  we  may 
predict  some  as  yet  unknown  future  point  in  the  series  from 
the  data  available.   A  very  simple  analogy  can  be  given.   If 
the  sales  in  a  certain  business  fluctuate  from  month  to 
month  as  shown  in  Figure  3*2.  and  we  compare  each  month's 
sales  with  the  preceding  month's  sales,  we  can  see  that  fol- 
lowing every  high  sales  month  there  is  a  lower  sales  month, 
and  vice  versa. 


-t- 
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Figure    3«2      Sales    Function 
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It  is  very  simple  to  predict  that  in  the  future  for  every 
high  month,  there  will  be  a  lower  month  to  follow.   However, 
if  we  increase  our  "lag"  and  compare  only  every  second  month, 
we  find  that  it  is  harder  to  make  such  predictions.   If  we 
go  further  and  compare  only  the  third  months,  we  lose  even 
more  confidence  in  our  predictions.   This  is  a  very 
simplified  explanation  of  Figure  3-3  which  generally  depicts 
the  autocorrelation  function. 


©     i 

Lag  in  months 


Figure  3*3  Autocorrelation  Function 


Figure  3*3  shows  that-  the  autocorrelation  function  is 
maximum  at  a  minimum  lag  and  generally  decreases  as  the  lag 
is  increased.   Relative  maxima  in  the  function  appear  which 
correspond  to  basic  cycles  existing  in  the  data,  but  unless 
the  function  is  periodic,  i.e.,  repeats  itself  identically 
at  some  latter  period,  the  value  of  the  autocorrelation 
function  will  never  approach  the  value  determined  ,at  zero 
lag.   Likewise,  as  the  lag  increases,  our  confidence  of 
prediction  decreases. 
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While  we  will  not  make  use  of  spectral  analysis  or  the 
autocorrelation  theory  to  separate  and  analyze  signals,  make 
predictions,  or  even  use  the  correlation  properties  of  the 
latter,  we  will  attempt  to  take  advantage  of  the  properties 
of  each  to  characterize  or  represent  functions  which  are 
generated  by  our  porous-media  samples.   For  example,  the 
autocorrelation  function  has  a  property  of  retaining  all  of 
the  harmonics  of  a  given  function,  containing  no  new  ones, 
and  discarding  the  phase  angle.   The  usefulness  of  this 
property  shall  become  apparent  as  we  explain  the  application 
of  these  statistical  methods  to  the  work  of  characterizing 
porous  media. 

With  the  background  definitions  of  the  power  spectrum 
and  the  autocorrelation  theory,  we  may  now  show  the  mathe- 
matical derivation  of  each. 

C .   Periodic  Functions 

52 
Lee   shows  that  a  periodic  function  may  be  expressed 

as  a  Fourier  series,  as  discussed  in  Appendix  A. 


f(x)  =  %  ♦   E   a„  cos  ^  +  b  sin  ^        (3-D 
2    n-1   n      L     n      L 


where: 


a  ■-  r  '   '  f(x)  c°s  ^  dx  (3.2 

n    L   -L/2  L 
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P   +L/2  nny 

b   -  f   S  f(x)  sin  Zy2-   dx  (3-3) 

.   n    L   -L/2  L 

L  ■  Fundamental  period  of  the  function. 

x  =  value  of  the  independent  variable  at  equally 

spaced  points  along  the  function,  such  that  x  goes 

from  -L/2  to  L/2. 
n  ■  number  of  harmonic  terms  to  which  the  series  is 

expanded . 
Equation  (3*1)  may  be  rewritten: 

n  =o°        .   x 
f(x)  -   L   F(n)e  ^L  (3.^) 

n=_oo 

in  which: 

F(n)  =  1  (a  -jb  );  for  n  -  0,  1,  2,  .  . .   (3-5) 
2 

j  =  the  imaginary  operator,  equal  to  the  square  root 
of  minus  (negative)  one.   (Some  texts  use  the 
notation  "i,"  and  a  more  definite  description  may 
be  found  in  any  text  dealing  with  complex  variables.) 

combining  (3«2)  and  (3*3) 

0/  \    1  r  +  L//2  w  v   -jnnf  dx  (3.6) 

F(n  )  -  t-  J  f  (x)  e     L 

L   -L/2 

1- 

Equation  (3«6)  is  the  Fourier  Transform  of  the  periodic 
function.   Lee  shows  that  this  equation  represents  an 
"analysis"  of  the  amplitudes  and  phase  angles  of  each  sinusoid 
into  which  function,  f('x),  is  resolved.   F(n)  is  also  called 
the  complex  spectrum  for  this  reason. 
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To  separate  the  amplitude  and  phase  characteristics  of 
the  complex  spectrum,  F(n),  we  may  write  it  as: 


(n)  =  i  JaJ  +  b2 


n   n 


e  xp 


j  tan 


-1 


n 


(3-7) 


where  the  amplitude  is  depicted  by   F(n)|  : 
and  the  term  d  ,  indicates  the  phase  angle: 


(3-8) 


On  ■  tan 


-1 


a 


(3-9 


Combining  the  above,  Equation  (3»^)  may  be  written: 


f(x) 


n' 


n  --00 


iX^ 


n   n 


cos 


nnx  4 
L     n 


(3.10) 


To  analyze  possible  correlation  between  two  functions, 
we  use  the  function,  t  -,  9  »  where: 


*«L/2 

T  n  o  cv  )  -  f  !        M*)  f9(x+r)  dx 
1,2  l    _L/2    1    l 


(3.11) 


We  call  the  left-hand  member  of  Eauation  (3.11)  the 
correlation  function  at  a  given  lag  ( T*  ) .   By  multiplying  one 
function  by  the  other,  at  various  values  of  Y  ,  we  can  deter- 
mine the  points  of  best  correlation. 

If  we  adopt  the  special  condition,  where  f-i(x)  "'  fp(x), 
we  can  write  (3.11)  as: 
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*  1,1  <*> 


I  / 


L/2 


fn  (x)  f.  (x  +  T  )  dx 
L   -L/2   X      1 


(3-12) 


Lee  shows  that  the  autocorrelation  function  is  also 
equat  to: 


n"-« 


F(n) 


2eKr 


(3-13) 


When  *Y  equals  zero,  Equation  (3-13)  reduces  to: 


*1.1  ^ 


,   +L/2 

—  /      f  ' 

L   -L/2   ! 


(x)dx 


CD 


n 


F1(n 


(  3  •  1^  ) 


or  the  mean  square  value  of  the  function  equals  the  sum,  for 
all  harmonics,  of  the  square  of  the  absolute  value  of  the 
spectrum.   The.  power  spectrum  is  then  defined  by: 


4>"  (n) 


F1(n) 


(3.15) 


It  can  be  demonstrated  that  the  power  spectrum  and  the 
autocorrelation  function  are  Fourier  transforms  of  each  other. 
As  previously  stated,  all  functions  consisting  of  identical 
sinusoidal  components  have  the  same  autocorrelation  function, 
regardless  of  phase  angle.   Thus,  the  power  spectrum,  as  a 
Fourier  transform  of  the  autocorrelation  function,  is  like- 
wise independent   of  initial  phase  angle  and  having  determined 
the  latter,  the  former  is  uniquely  determined. 

Lee  also  shows  that  these  functions  may  be  written  in 
terms  of  Fourier  coefficients: 
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If  |Pj.(n)|-  |(an'Jbn) 
then : 


<h  ,    v  .  |«/  x|  2  .a  2+b  2  (3.16) 

vp  (n)  ■   F(n)    "  n    n  ** 

1    '    5 

and: 

*  lfl  (T  >  -  !°!  +  I  £  (V  V)  c°s  T^  (3.17) 

From  Equations  (3.16)  and  (3-17) »  it  is  possible  to 
calculate  these  desired  properties  of  a  function  with  the  use 
of  a  computer  program.   (See  Appendix  C.) 

D.  Aperiodic  Functions 

Whereas  the  tool  of  analysis  for  handling  a  periodic 
function  is  the  Fourier  expansion,  it  is  necessary  to  rely 
upon  the  Fourier  integral  to  analyze  the  aperiodic  function, 
The  derivation  of  this  integral  and  the  steps  necessary  to 
evaluate  the  power  spectrum  and  autocorrelation  function  for 
this  type  of  function  are  not  required  within  the  scope  of 

this  paper.   A  discussion  of  the  theory  is  contained  in  texts, 

53        5*4- 
such  as  Tolstov^  or  Lee. 

E.  Random  Functions  » 

We  define  a  random  function  as  follows: 

1.   The  generating  process  cannot  be  completely  evaluated 
and  a  precise  and  unique  determination  of  the  value 
of  the  dependent  variable  for  any  value  of  the 
independent  variable  cannot  be  made. 
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2.   As  a  conseauence  of  (1),  the  value  of  the  dependent 
variable  cannot  be  predicted  at  any  given  future 
point  in  time. 
In  other  words,  a  random  function  is  one  that  is  generated 
by  some  unknown  source  which  dictates  the  manner  in  which  the 
function  behaves,  but,  due  to  its  stochastic  or  probabilistic 
nature,  cannot  provide  the  actual  value  of  the  function  at 
any  given  moment  in  time.   An  example  of  a  random  function 
might  be  the  noise  that  a  radio  picks  up  in  addition  to  the 
intended  message  or  program.   This  noise  may  be  the  results 
of  thermal  emissions,  static,  or  spillover  from  other  fre- 
quencies, but  an  exact  description  of  the  noise  is  not  pos- 
sible because  of  the  complexity  and  random  nature  of  the 
noise  source.   It  might  be  noted  that  in  most  economic 
series,  although  containing  obvious  trends  caused  by  seasonal 
or  secular  variations  which  may  be  isolated,  certain  strong 
random  components  exist. 

Because  of  the  problems  created  by  this  "noise,"  it  is 
of  great  import  in  communications  work  that  analysis  of  such 
functions'  be  a  possibility.   In  order  to  analyze  random 
functions,  we  find  it  necessary  to  introduce  certain  re- 
strictions, and  assumptions.   First,  we  assume  that  the 
function  has . an  infinite  period  and  that  during  this  period, 
the  function  will  trace  a  pattern,  such  that  the  value  of 
the  dependent  variable  will  pass  through,  or  infinitely 
close  to,  all  possible  values  which  are  consistent  with 
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the  "ergotic  hypothesis"  as  espoused  by  Maxwell,  and 
discussed  in  Lee.   Accepting  this  theorem,  we  may  consider 
our  random  generating  process  as  generating  a  group  of 
finite  length  random  functions  rather  than  one  infinitely 
long  function  having  the  combined  characteristics  of  the 
entire  grouping.   The  individual  functions  follow  the  nature 
of  the  infinite  function  on  an  average  value  basis  and  in  a 
probabilistic  manner.   By  the  latter,  it  is  meant  that  such 
properties  as  mean  value  crossings,  or  zero  value  crossings 
per  unit  length  follow  the  probability  distribution  of  the 
parent  function.   These  finite  functions  will  be  referred  to 
as  an  ensemble  of  random  functions.   They  are  all  generated 
by  the  same  source.   This  grouping  of  the  ensemble  is  a  very 
necessary  step  in  the  analysis  of  random  functions,  as 
realized  and  set  forth  by  Gibbs .  ^ 

In  such  an  ensemble  of  functions,  each  of  the  finite 
functions  will  .be  called  member  functions.   It  is  important 
to  realize  that  we  assume  that  each  of  these  member  functions 
is  generated  by  the  same  source,  although  we  cannot  mathe- 
matically formulate  the  operation  of  the  source.   This  being 
so,  we  sense  that  properties  which  characterize  any  member 
of  the  ensemble,  in  turn  characterize  all  others,  as  well  as 
the  infinite  function  from  which  the  ensemble  was  created. 
One  of  these  properties  would  be  the  autocorrelation  function, 

Hn,  +L 

^      (T)  -  lim    L  /   f1(x)fn(x*f)dx    (3.18) 
1 ' l  L— *  »   2    -L  l  1 
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which  we  assume  exists  for  all  values  of T'  .   Another  such 
property  is  the  power  spectrum.   Wiener  '  developed  the 
theorem  that  the  autocorrelation  function  of  a  random  func- 
tion and  the  power  density  spectrum  are  related  to  each 
other  as  Fourier  Cosine  transforms.   The  fundamental  assum- 
ption underlying  this  theorem  is  that-  the  autocorrelation 

function  exists  and  is  common  for  all  member  functions  of 

57 
the  ensemble.   Blackman  and  Tukey   state  that,  although 

Equation  (3. 18)  is  frequently  called  the  autocorrelation 

function,  a  more  appropriate  name  for  the  term  defined  by 

Equation  (3«18)  is  the  autocovariance  function.   They  state 

that  the  autocorrelation  function  should  be  defined  as  the 

normalized  autocovariance  function  obtained  as  the  ratio 

of  the  latter  function  at  any  lag  to  the  function  at  zero 

lag.   This  affords  the  autocorrelation  function  the  property 

of  having  the  value  of  unity  at  zero  lag.   However,  it  is 

convenient  for  our  purposes,  and  entirely  consistent  with 

the  definitions,  to  use  the  names  interchangeably  for 

characterizing  porous  media. 

In  order  that  convenient  equations  might  be  developed 

which  will  allow  us  to  calculate  these  properties  with  the 

•  ■  •  CO 

aid  of  a  computer,  we  refer  to  Granger  and  Hatanaka.  They 

show  that  for. real  processes  the  autocovariance, >y~  ,  may  be 
written: 

Hy       ■  2/^  cos  r  uu  f(iD)dU)  (3.19) 

where:      uu  ■  =%•     <   o 
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Formal  inversion  provides  the  power  spectrum 

f(u>)  -  iy,  (u/Q    +   1  V.  cos  juuj  (3.20) 

For  finite  data,  using  numerical  integration,  an 
approximation  to  Equation  (3.19)  is: 

•r-'.k.  ^(^i-*)(^r-)  (3'21). 

where  n  =  number  of  data  points 

The  estimate  of  Equation  (3.20)  becomes 

f(u,)  =  2w(C°  *  2   j^  °0  °OS  ^)         •'      (3.22) 

Hannan  shows  that  this  is  a  poor  estimate,  but  that 
it  may  be  improved  by  properly  weighting  the  values  of  the 
autocovariance  terms,  so  that 

•a  /  m-i  \ 

f  (u))  =  ^-  I   cQXo(uu)  +  2   2   c^A^uu)  cos  jwj     (3-23) 

where  "m"  is  arbitrarily  selected  by  the  user  and  de- 
notes the  number  of  lags  used.   It  should  be  less  than 
one-third  the  value  of  "n . "   X.(w)  represents  the  value  of 

the  weighting  factors. 

There  is  considerable  discussion  concerning  the  choice  • 
of  the  weighting  factors,  X(uu).   However,  the  two  most 
widely  used,  developed  by  J .  W.  Hamming  and  J.  Von  Hann , 
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are  very  similar  and  for  this  paper  we  will  use  the  latter's 
system,  now  commonly  known  as  the  Tukey-Hanning  estimate. 
It  has  the  values 


*  -iti'-F] 


(3.2*0 


The  actual  equations  for  estimating  the  covariance  and 
the  power  spectrum  may  now  be  written 
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(3.26) 


U.  =.25LH  +  .5V  .25L.+1 


(3.27) 


where  L  ,  +  L.,  ,  L  .,  4  L  -, 
-1    41'   m*l    m-1. 


The  L.  terms  are  the  raw  estimate  of  the  spectrum,  and 

J 

the  U.    terms  are  the  smoothed  values.   The  smoothing 
operation  of  Eauation  (3.27)  is  apparent  from  the  form  of 
the  equation. 

To  illustrate  Equations  (3.25)  and  (3.26),  we  might 
examine  two  very  simple  functions.   Function  A  is  an  on-off 
function  as  shown  in  Figure  3«^. 
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Figure  3«^   Example  Function  A 


By  examining  the  figure,  it  is  seen  that  the  function 
repeats  every  eight  units  and  that  the  autocovariance 
function  has  the  shape  shown  in  Figure  3»5« 


+  1 


-  I 


8    »o     12.  14- 


it, 


Figure  3»5   Autocovariance  Estimate  of  Function  A 


Referring  to  Equation  (3.26),  it  can  be  seen  that  the 
power  spectrum  will  have  a  peak  value  when  the  frequency 


parameter  "j"  is  such  that  the  cosine  term,  cos 


m 


,  has  a 


period  coincident  with  the  period  of  the  autocovariance 
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function.   As  an  example,  taking  m  -  150  lags,  the  power 
spectrum  reaches  its  peak  value  when  the  frequency  parameter 
" j"  is  equal  to  37 • 5 •   The  power  spectrum  appears  as  shown 
in  Figure  3*6. 


15 
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Figure  J. 6      Power  Spectrum  of  Function  A  at  150  Lags 


These  results  are  easily  checked  by  observing  that  the 
original  function  passes  through  37'"   cycles  in  150  units. 

A  second  function  which  can  be  used  to  demonstrate  the 
nature  of  the  analysis  performed  by  use  of  the  auto- 
correlation function  and  the  power  spectrum  is  Function  B. 
Function  B  is  a  combination  of  two  simple  on-off  functions 
which  have  been  superimposed.   The  frequency  of  the  first 
function  is  four  cycles  per  eighty  points.   The  frequency  of 
the  second  function  is  ten  cycles  per  eighty  points. 
Function  B,  as  shown  in  Figure  3«7»  has  very  little  apparent 
regularity. 
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Figure  3*7   Example  Function  B 


The  power  spectrum  for  50  lags  shows  peaks  at  5  and 
12. 5>  which  were  the  individual  peaks  for  each  basic  fre- 
quency used  in  the  combined  function.   In  other  words, 
Equation  (3.26)  effectively  selects  the  fundamental  fre- 
quencies which  constitute  a  function,  even  though  they  may 
not  be  apparent  individually.   Further,  a  measure  of  the 
relative  contribution  of  each  frequency  is  offered. 

The  computer  program  "COVAR,"  discussed  in  Appendix  C, 
has  been  written  specifically  to  perform  the  mathematical 
operations  set  forth  in  Equations  (3*25)  and  (3*26) .   The 
above  figures  were  made  from  data  and  results  of  that  pro- 
gram, run  on  Functions  A  and  B. 


CHAPTER  IV 
EXPERIMENTAL  PROCEDURE 

A.  General 

In  /this  thesis  an  attempt  is  being  made  to  classify 
or  characterize  porous  media.   It  is  our  goal  to  differ- 
entiate between  media  which  are  not  statistically  identical 
through  the  use  of  various  forms  of  statistical  analysis. 
In  order  to  apply  the  various  tools  outlined  in  Chapter  III, 
we  require  a  variety  of  samples  which  we  know  to  be  disimilar. 
In  this  chapter,  we  will  discuss  the  selection  and  preparation 
of  these  samples  and  the  collection  of  data  therefrom. 

B.  Preparation  of  Samples 

Data  used,  in  this  thesis  were  taken  from  two  types  of 
porous-media  samples.   The  first  sample  was  a  photomicrograph 
of  an  actual  porous  medium.   The  other  samples  used  were 
prepared  by  the  author  by  imbedding  mult i -size  spherical 
glass  beads  in  plastic. 

The  photomicrograph  is  a  photograph  of  an  actual  porous 
medium  taken  through  a  microscope.   In  the  sample  used,  the 
photographic  image  wa§  enlarged  ninety-two  times  from  an 
original  size  of  approximately  .12  centimeters  by  .16 

centimeters.   This  photograph  was  taken  from  the  paper  by 
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Fara  and  Scheidegger.    The  original  medium  was  sandstone 

with  a  porosity  of  26.3  percent.  % 

3^ 


3  5 
The  other,  samples  used  were  created  by  the  author  in  the 
following  manner: 
Materials  : 

(1)  Laboratory  glass  beads,  three  sizes: 
3mm.  ,  4mm.,  5mm. 

(2)  Natco  Premium  Quality  Clear  Polyester  Casting 
Resin,  NL-600. 

(3)  Violet  dye  (water  base). 

(4)  Waxed  paper  cups. 
Procedure : 

(1)  Mix  dye  and  resin  to  achieve  a  colored, 
slightly  opaque  plastic. 

(2)  Place  liquid  plastic  in  waxed  paper  container. 

(3)  Drop  glass  beads,  in  varying  proportions  as 
required,  into  the  plastic. 

(4)  Place  'in  over  for  one  day  (30  degrees  centi- 
grade).  Allow  to  set. 

(5)  Cut  through  sample  in  a  random  manner  to 
realize  an  unbiased  cross  sectional  surface 
of  the  mold.   This  surface  area  represents  a 
model  of  a  porous  medium.   Where  the  glass 
beads  are  apparent  is  taken  as  solid space, 
and  where  plastic  is  apparent  is  taken  to  be 
pore  space. 


36 
For  this  paper,  three  samples  were  prepared.   Sample  X 
consisted  of  3mra .  beads  only.   Sample  Y  was  a  mixture  of 
equal  volumes  of  4mm.  and  5mm*  beads.   Sample  Z  was  a  mixture 
of  equal  volumes  of  3mm.,  4mm.,  and  5mm.  beads.   The  cut 
samples  were  approximately  two  inches  in  diameter. 

C .   Preparation  of  Data 

Following  the  preparation  of  the  samples,  photographs 
of  each  were  made  by  the  Photo  and  Graphic  Arts  Department, 
University  of  Kansas  (see  Appendix  D).   The  images  of  the 
created  samples  were  expanded  six  times,  while  that  of  the 
Scheidegger  photomicrograph  was  already  of  appropriate 
dimensions.   An  80  x  80  grid,  with  ten  units  per  inch  was 
overlaid  on  each  photograph.   At  each  grid  point,  the  data 
in  point  was  recorded  as  being  a  41  if  a  pore  space  (plastic), 
and  a  -1  if  a  solid  space  (beads).   These  80  x  80  data 
points  were  easily  converted  to  input  information  on  eighty 
column  computer  cards,  for  use  with  the  programs  written  for 
this  thesis.   (The  negative  values  were  read  into  the 
computer  as  fixed  point  zeros  and  changed  internally  to  their 
correct  values  as  negative  ones.) 


CHAPTER  V 
DISCUSSION  OF  RESULTS 

The  mission  of  this  project  was  to  provide  means  of 
statistically  analyzing  porous-media,  to  show  statistical 
variation  in  media  known  to  be  different,  and  to  offer 
direction  and  inspiration  to  further  work  in  this  area.   The 
author  fully  realizes  that  from  the  limited  number  of  sam- 
ples dealt  with  in  this  thesis,  there  is  insufficient  evidence 
to  prove  that  the  far  reaching  goal  of  statistically  charac- 
terizing all  media  can  be  accomplished.   However,  in  this 
chapter,  we  will  offer  indications  that "the  above  outlined 
missions  have  been  successfully  completed  and  that  success 
of  the  overall  goal  mentioned  is  entirely  possible.   It  is 
hoped  that  from  the  results  shown,  others  will  be  stimulated 
to  continue  research  toward  that  end. 

Two  basic  computational  tools  were  used  to  statisti- 
cally analyze  the  data  taken  from  the  porous-media  samples. 
Although  discussed  in  Chapter  IV,  it  might  be  well  to  point 
out  again  that  these  data  were  taken  at  equally  spaced  points 
along  lines  randomly  drawn  across  the  cut-face  of  the  samples. 
Where  the  line  passed  thrbugh  pore  space,  a  value  of  41  was 
assigned  to  the  function,  and  where  the  line  bisected  solid, 
a  -1  value  was  assigned.   A  series  of  these  points  taken  in 
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order  along  the  line  then  formed  a  simple  orthogonal  function 
with  the  +1  to  -1  cycle  arranged  in  a  seemingly  random  pat- 
tern.  Based  upon  the  assumption  that  porous  media  are 
statistically  constituted,  it  was  hoped  that  several  of  these 
lines  or  functions  taken  from  the  same  medium  would  be 
statistically  similar  while  those  taken  from  different  media 
would  vary  by  what  is  known  as  a  "significant  difference." 

The  first  method  used  to  statistically  analyze  the 
functions  was  to  fit  each  with  a  single  Fourier  series  using 
"ID  Fourier"  computer  program.   (See  Appendixes  A  and  C.) 
Although  each  function  was  considered  as  being  random  in 
nature,  it  was  analyzed  over  a  finite  length.   The  Fourier 

coefficients  were  determined  and  then  combined  according  to 

2     2     2 
the  Lee  estimate  of  the  power  spectrum,  c    -  a   +  b   ,  so 

f      ^      '   n     n     n  ' 

that  the  initial  starting  point  or  phase  angle  might  be  dis- 
regarded.  These  values  of  c  ,  as  estimates  of  the  power 
spectrum,  are  measurements  of  the  contributions  of  the 
fundamentals  frequencies  constituting  each  function.   This 
is  essentially  the  technique  used  in  the  paper  by  Fara  and 
Scheidegger. 

Typical  results  of  the  single  Fourier  series  approach 
are  displayed  in  Table  5*1  Shown  are  the  c  values  for  four 
lines  across  the  same  media.   As  can  be  readily  observed, 
there  is  very  little  consistency  among  the  values,  indi- 
cating that  this  method  yields  a  poor  estimate  of  the  power 
spectrum.   As  a  result,  the  simple  Fourier  approach  was  not 
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pursued  to  any  great  detail,  although  later  in  this  chapter 
we  will  compare  the  results  of  using  this  approach  with  those 
of  our  second  technique  on  an  identical  sample. 

Because  of  the  poor  results  achieved  in  the  single 
Fourier  approach,  no  analysis  of  the  two  dimension  cross 
sectional  faces  of  the  media  was  attempted.   However,  a 
computer  program  "2D  FOURIER"  to  determine  double  Fourier 
series  coefficients  for  a  surface  was  developed,  tested,  and 
is  fully  operational.   (See  Appendixes  B  and  C.) 

The  second  statistical  tool  used  to  analyze  the  data 
was  the  determination  of  the  communication  theory  properties 
of  the  functions  or  lines.   These  properties,  the  autoco- 
variance  function  and  the  power  spectrum,  are  discussed  in 
Chapter  III.   The  properties  were  calculated  using  Fortran 
equivalents  of  Equations  (3.25)  and  (3.26)  in  the  computer 
program  "COVAR."   (See  Appendix  C.)   As  a  part  of  that  pro- 
gram, the  power  spectrum  was  plotted  for  each  of  several 
lines  per  media  type,  as  well  as  the  90  and  95  percent 
confidence  bands.   These  bands  were  created  by  multiplying 
the  smoothed  power  spectrum  by  constants  presented  in  Granger 
and  Hatanaka.    The  constant  used  is  a  function  of  the 
"equivalent  degrees  of  freedom,"  an  approximation  chosen  to 
allow  for  the  fact  that  the  series  under  consideration  is 
not  necessarily  exactly  normally  distributed.   These  bands 
indicate  90  or  95  percent  confidence  that  the  true  value  of 
the  power  spectrum  lies  within  that  range  of  values.   In 
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other  words,  90  or  95  percent,  as  the  case  may  be,  of  all 
power  spectra  values  from  the  same  population  should  lie 
within  that  range.   Generally,  the  approach  was  to  compare 
data  lines  from  the  same  media  to  see  if  the  observed  plots 
of  the  power  spectra  were  similar.   Then  it  was  determined 
whether  the  values  of  each  function  within  the  sample,  at 
each  frequency,  fell  within  the  confidence  bands  of  the 
other  lines  within  that  sample  or  type  of  media.   A  similar 
technique  was  used  to  compare  the  power  spectra  profile  of 
one  sample  media  with  that  of  another  media. 

For  these  comparisons,  the  effect  of  sample  lengths  was 
investigated.   Samples  of  various  lengths  were  created  by 
combining  the  eighty  point  crossings  developed  in  accordance 
with  Chapter  IV.   The  crossings  were  shuffled  to  reduce 
bias  and  eliminate  peaks  which  would  occur  in  the  autoco- 
variance  function  at  eighty  point  intervals.   These  peaks 
would  be  due  to  the  similarity  between  adjacent  lines  drawn 
across  a  porous  media.   It  was  discovered  after  a  series  of 
trials  that  the  longer  the  lines,  the  less  bias  there  was 
to  the  idiosyncrasies  of  single  crossings,  i.e.,  a  line 
might  pass  through  continuous  solid  or  pore  space  for  an 
many  as  twenty  or  thirty  points  without  interruption  or 
change.   Longer  lines  tend  to  minimize  the  effect  of  these 
rare  exceptions.   As  a  result,  lines  of  3200  data  points  or 
^0  crossings  were  utilized,  this  being  the  maximum  number 
of  points  that  could  be  handled  on  the  computer  at  one  time. 
It  might  be  noted  that  lines  of  1600  data  points  afforded 
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very  similar  results,  indicating,  that  the  lines  of  3200  points 
were  of  sufficient  length. 

It  was  discovered  that  the  number  of  lags  selected,  or 
more  properly  stated,  the  length  of  the  maximum  lag  selected, 
played  a  very  important  part  in  the  analysis.   Again  after 
many  trials,  it  became  apparent  that  when  the  maximum  lag 
selected'  was  too  short,  the  apparent  frequencies  were 
"packed"  at  the  very  low  power  spectrum  values.   In  other 
words,  with  many  frequencies  constituting  the  function, 
each  band  or  grouping  of  frequencies  was  too  wide  or  contained 
too  many  frequencies  to  effectively  delineate  the  function. 
To  long  a  maximum  lag  had  the  effect  of  causing  random  spikes 
in.  the  power  spectrum  which  were  characteristic  of  an  in- 
dividual line 'but  not  of  the  media  as  a- whole.   These  "spikes" 
might  be  considered  a  result  of  sampling  variability.   Also, 
at  long  lags,  values  of  the  autocovariance  function  tended 
to  distort  the  .power  spectrum  profile.   In  light  of  the. 
above,  a  lag  of  40  was  selected  as  being  optimum.   It  had 
the  effect  of  smoothing  the  data  and  still  offering  good 
disorimination . 

Figures  5*1  through  5«3  offer  three  plots  of  the  power 
spectra  of  different  media.   Appendix  E  contains  the  com- 
puter results  from  which  these  lines  were  plotted  and  the 
computer  plots.   The  three  profiles  are  obviously  different, 
but  it  is  important  to  discuss  these  differences  in  a  more 
absolute  statistical  manner.   As  shown  in  the  plots,  the 
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95  percent  confidence  band  of  one  of  the  lines  is  super- 
imposed upon  the  figure.   By  looking  at  each  figure,  it  can 
be  seen  that  generally  the  power  spectrum  of  the  other,  or 
"outer"  line  lies  within  that  confidence  band.   More 
specifically,  the  outer  line  is  within  the  band  at  24  points 
out  of  40  in  Figure  5«1«   This  is  the  equivalent  of  60  per 
cent  of  the  time.   In  Figure  5.2  the  result  is  23  out  of  40 
(57-5%)   and  in  5*3  it  is  again  23  out  of  40  {57-5%). 

While  these  results  violate  the  95  percent  criteria., 
it  is  important  to  note  that  many  of  the  points  are  very 
close  to  the  limits  and  as  only  two  examples  were  displayed 
per  sample,  variability  of  samples  has  a  great  effect. 
Ideally,  many  more  lines  per  sample  should  be  obtained, 
plotted  and  analyzed.   Then  the  confidence  bands  around  the 
mean  of  those  lines  should  be  superimposed  upon  a  plot  of 
the  lines.   Looking  at  all  of  the  power  spectra  values  that 
lie  within  these  bands,  the  95  per  cent  criteria  would  be 
approached.   Using  this  method  on  four  lines  of  shorter 
length  from  these  same  media,  results  very  close  to  that 
criteria  were  realized  and  as  more  samples  are  considered, 
it  is  expected  that  the  95  percent  criteria  would  be  ful- 
filled.  However,  this  was  not  proven  in  this  work.   By  the 
time  the  tools  to  analyze  media  had  been  developed,  and  made 
available  in  the  form  of  fully  tested  computer  applications, 
and  when  the  results  discussed  above  became  apparent,  there 
was  insufficient  time  to  create  or  acquire  additional  samples 
for  this  thesis. 


43 
Cross  comparing  samples,  there  is  very  little  correla- 
tion in  comparison  to  the  correlation  within  the  same  media. 
Between  Figures  5»1  and  5»3»  which  are  very  similar  media 
(see  Chapter  IV),  the  base  lines  are  within  the  other's 
confidence  bands  only  12  of  40  points,  or  30  percent  of 
the  times.   Between  Figures  5.1  and  5«2,  the  result  is  only 
one  out  of  40  ( 2 . 5% )  ,  and  between  5*2  and  5 • 3 »  it  is  one  out 
of  k§    [2.5%),      Visually  comparing  the  different  spectra,  it 
can  be  seen  that  the  3mm .  beads  sample,  which  only  has  small 
beads  in  its  make-up,  has  lower  values  than  the  other 
samples  at  low  frequencies,  and  higher  values  than  the  others 
at  higher  frequencies.   This  is  to  be  expected  as  there  is 
a  greater  fluctuation  in  the  41  to  -1  cycle  in  the  original 
3  mm.  bead  data.   There  is  less  opportunity  for  long  strings 
of  41  or  -1  values  due  to  this  composition  of  the  sample. 
Similar  analysis  can  be  made  between  the  other  samples. 

It  is  felt  then,  that  this  is  a  definite  indication  that 
the  general  shape  of  the  power  spectrum  shown  is  indicative 
of  the  composition  of  the  media  from  which  the  data  is 
taken.   Spectra  from  different  media  are  not  only  visibly 
different,  but  also  it  is  noted  that  the  differences  appear 
to  be  statistically  significant.   Again  it  is  repeated  that 
many  more  samples  of  various  types  need  to  be  analyzed  to 
prove  this  conclusively,  but  the  technique  is  readily  avail- 
able, and  indications  are  that  the  results  would  be 
favorable. 


44 
Of  interest,  Table  5*2  and  Figure  5-^  display  the  com- 
puted autocovariance  estimate  and  the  power  spectrum  of  the 
actual  rock  sample  discussed  in  Chapter  IV. 

A  brief  comparison  was  made  between  the  two  techniques 
used  to  analyze  the  functions,  the  single  Fourier  series  and 
the  power  spectrum.   Figure  5*5  shows  this  comparison.   The 
function  analyzed  was  the  line  drawn  across  the  porous  media 
shown  in  the  Fara-Scheidegger  paper.   (See  Chapter  IV).   As 
can  readily  be  observed,  there  is  very  little  consistency 
in  the  results.   This  is  further  indication  that  the  single 
Fourier  series  approach  is  not  particularly  suited  for  this 
type  of  analysis. 
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FIGURE  5.*( 


POWfR  SPhCIRUM  ESTIMATION 

FOR  DISCRfTF  TIME  SERIES 

USING  TUKEY  MANNING  WEIGHTING  FACTORS 


SCHEIOEGGER  DATA   1?00  POINTS 


V' 


SMOOTHEO  POWER  SPECTRUM 
WITH  90  PCT.  CONFIDENCE  LIMITS 
I     2  3     *     b  6     T 


SMOOTHEO  POWER  SPECTRUM 
WITH  95  PCT.  CONFIDENCE  LIMITS 
12     1*967 
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TABLE  5.1 

Tabulated  Results  of  Single  Fourier  Series 
Analysis  of  Three  Size  Beads  Data 


Power  Spectrum. 

Estimate 

1 

2 

3 

4 

C(0) 

.2925 

.5637 

.4425  ' 

.2900 

C(l) 

.0005 

.0012 

.0048 

.0298 

C(2) 

.0033 

.0002 

.0039 

.0029 

C(3) 

.0057 

.0000 

.0001 

.0029 

C(4) 

.0020 

.0002 

.0007 

.0006 

c(5) 

.0078 

.0002 

.0017 

.0023 

C(6) 

.0002 

.0002 

.0017 

.0023 

C(?) 

.0028 

.0002 

.0023 

.0140 

C(8) 

.0058 

.0000 

.0041 

.0046 

C(9) 

.0022 

.0109 

.0018 

.0036 

C(10) 

.0027 

.0028 

.0017 

.0005 

C(ll) 

.0020 

.0011 

.0151 

.0071 

C(12) 

.0012 

.0016 

.0008 

.0018 

C(13) 

.0030 

.0000 

.0080 

.0050 

C(l4) 

.0005 

.0006 

.0040 

.0001 

C(15) 

.0034 

.0003 

.0049 

.0018 

C(l6) 

.0049 

.0000 

.0012 

" .0024 

■C(17) 

.0048 

.0007 

.0017 

.0079 

C(18) 

.0049 

.0003 

.0007 

.0017 

C(19) 

.0202 

.0006 

.0087 

.0187 

C(20) 

.0003 

.0120 

.0048 

.0056 

C(30) 

.0024 

.0124 

*  .0012 

.0015 

C(40) 

.0047 

.0139 

.0076 

.0039 

c(50) 

.0119 

.0066 

.0148 

.0010 

C(0)  is  an  indication  of  the  porosity  of  the  media.   The 
actual  porosity  of  the  sample  was  35 • 65  percent. 
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POWTR  SPlCfRUK  ESTIMATION 
TABLE  5.2  POM  OISCREIF  TIMF  SERIFS 
USIKG  TUKFV  MANNING  WEIGHTING  FACTORS 


SCHEIOEGGER  DATA   3200  POINTS 


COVARIASCE 

RAM 

LAG 

SPECTRUM 

SPECTRUM 

0 

O.OOGO 

0.653B 

1 

0.00*00 

0.7803 

2 

0.0000 

0.6303 

3 

o.ooco 

0.4465 

A 

0.0000 

0.5051 

5 

0.0000 

0.42*7 

6 

040000 

0.3666 

7 

0.0000 

0.2889 

8 

0.0000 

0.2249 

9 

0.0000 

0.1705 

10 

0.0000   — 

—  0.1755 

11 

0.0000 

0. 1 170 

12 

o.ooco 

0. 1201 

I  J 

0.0000 

0.0962 

1* 

o.ooco 

0.0729 

15 

0.0000 

0.0857 

16 

o.ooco 

0.0614 

17 

0.0000 

0.0725 

IB 

0.0000 

0.0611 

19 

0.0000 

0.0556 

20 

0.0000 

0.0484 

?1 

0.0000 

0.0394 

71 

0.0000 

0.0419 

23 

0.0000 

0.0448 

2* 

0.0000 

0.0323 

25 

0.0000 

0.0336 

26 

0.0000 

0.0362 

27 

0.0000 

0.0327 

28 

0.0000 

0.0282 

29 

0.0000 

0.0235 

30 

0.0000 

0.0251 

31 

0.0000 

0.0279 

32 

0.0000 

0.0267 

33 

0.0000 

0.0273 

3* 

0.0000 

0.0262 

35 

0.0000 

0.0248 

36 

0.0000 

0.0193 

37 

0.0000 

0.0213 

38  . 

0.0000 

0.0286 

39 

0.0000 

0.0296 

*0 

0.0000 

0.0163 

SMOOTHED  SPECTRUM  WITH  INOICATEO  CONF 

90  PCT.  90  POT.  95  PCT. 

LOWFR                     UPPER  LOWER 

LIMIT  AVERAGE       LIMIT  LIMIT 

0.6482  0.7171  0.7716  0.6210 

0.6429  0.7112  0.7652  0.6159 

0.5622  0.6219  0»669l  0.5385 

0.45H4  0.5071  0.5457  0.4392 

0.4252  0.4704  0.5061  0.4073 

0.3890  0.4303  0.4630  0.3726 

0.32fO  0.1617  0.3892  0.3132 

0.2643  0.2923  0.3146  0.2532 

0.2055  0.2273  0.2446  0.1968 

n.1676  0.1853  0.1994  0.1605 

0.1443  0.1596  0.1718  0.1382 

0.1215  0.1344  J.1446  0.1164 

0.1061  0.1174  0.1263  0.1016 

0.0B89  0.0984  0.1058  0.0852 

0.0741  0.0819  0.0882  0.0710 

0.0691  0.0764  0.0823  0.0662 

0.0635  0.0703  0.0756  0.0608 

0.0604  0.0669  0.0719  0.0579 

0.0565  0.0625  0.0673  0.0542 

0.0449  0.0551  0.0593  0.0478 

0.0433  0.0479  0.0516  0.0415 

0.0382  0.0423  0.0455  0.0366 

0.0180  0.0420  0.0452  0.0364 

0.0370  0.0410  0.0441  0.0355 

0.0323  0.0358  0.0385  0.0310 

0.0107  0.0339  0.0365  0.0294 

0.0313  0.0347  0.0373  0.0300 

0.0293  0.0325  0.0349  0.0281 

C.0255  0.0282  0.0303  0.0244 

0.0227  0.0251  0.0270  0.0217 

0.0230  0.0254  0.0273  0.0220 

0.0243  0.0269  0.0289  0.0233 

0.0245  0.0271  0.0292  0.0235 

0.0243  0.0268  0.0289  0.0232 

0.0236  0.0261  0.0281  0.0226 

0.0215  0.0238  0.0256  0.0206 

0.0192  0.0212  0.0228  0.0184 

0.0205  0.0226  0.024*  0.0196 

0.0244  0.0270  0.0291  0.0234 

0.0235  0.0260  0.0280  0.0226 

0.0208  0.0230  0.0247  0.0199 


I0ENCE  LIMITS 


AVERAGE 

0.7171 

0.7112 

0.6219 

0.5071 

0.4704 

0.4303 

0.3617 

0.2923 

0.2273 

0.1853 

0.1596 

0.1344 

0.1174 

0.0984 

0.0819 

0.0764 

0.0703 

0.0669 

0.0625 

0.0551 

0.0479 

0.0423 

0.0420 

0.0410 

0.0358 

0.0339 

0.0347 

0.0325 

0.0282 

0.0251 

0.0254 

0.0269 

0.0271 

0.0268 

0.0261 

0.0238 

0.0212 

0.0226 

0.0270 

0.0260 

0.0230 


95  PCT. 
UPPER 
LIMIT 
0.7773 
0.7709 
0.6741 
0.5497 
0.5099 
0.4664 
0.3921 
0.3169 
0.2464 
0.2009 
0.1730 
0.1457 
0.1272 
0.1066 
0.0888 
0.0829 
0.0762 
0.0725 
0.0678 
0.0598 
0.0520 
0.0458 
0.0455 
0.0444 
0.0388 
0.0368 
0.0376 
0.0352 
0.0305 
0.0272 
0.0275 
0.0291 
0.0294 
0.0291 
0.0283 
0.0258 
0.0230 
0.0246 
0.0293 
0.0282 
0.0249 


CHAPTER  VT 
CONCLUSIONS  AND  RECOMMENDATIONS 

In  this  thesis,\an  attempt  has  been  made  to  provide 
effective  methods  of  statistically  analyzing  porous  media. 
It  has  been  shown  that  the  properties  derived  from  communi 
cation  theory,  the  autocovariance  estimate  and  the  power 
spectrum,  are  capable  of  differentiating  between  different 
media.  Further.,  it  would  appear  from  the  work  with  the 
simple  representations  of  media  which  were  created  for  this 
problem,  that  real  porous  media  lends  itself  to  this  type 
of  analysis. 

Perhaps  the  most  significant  conclusion  that  can  be 
drawn  from  this  thesis  then,  is  that  the  techniques  of 
statistical  communication  theory  are  adaptable  for  this  form 

of  analysis.  Although  the  classical  Fourier  approach  did 

63 

not  provide  consistent  results,  the  Tukey  Hanning  **  power 

spectrum  estimate  offered  good  discrimination  between  the 
different  media. 

There  are  several  obvious  steps  required  to  further  the 
proof  and  development  of  the  theory  presented  herein.  A 
large  number  of  samples  from  each  simple  type  of  media  should 
be  analyzed  to  affirm  the  consistency  of  results.  Then  new 
samples  should  be  fabricated  with  varying  proportions  of  the 
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53 
different  size  beads,  to  determine  how  fine  a  discrimination 
can  be  made. 

If  these  tools  continue  to  offer  good  results  with 
idealized  media,  the  next  step  would  be  to  test  them  on  real 
rock  samples.   Such  statistical  methods  as  Discrimination 
Function  Analysis  or  Factor  Analysis  might  be  applied  to  the 
data  to  improve  the  differentiation  between  "different" 
media. 

Finally,  if  these  techniques  prove  successful  upon  real 
media,  it  would  be  logical  to  then  attempt  to  relate  the 
statistical  properties,  such  as  the  power  spectrum,  to  the 
static  and  dynamic  fluid  behavior  of  the  med.ia. 
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APPENDIX  A 
SINGLE  FOURIER  SERIES 

In  many  technical  problems,  it  is  necessary  to  deal  with 
complicated  periodic  functions,  that  is  functions  F(x)  = 
F(x  +  T)  for  all  values  of  x,  T  being  th'e  period  of  the  func- 
tion.  Use  of  the  Fourier  series  allows  us  to  write  these 
complex  functions  in  terms  of  simple  trignometric  functions. 
Although  Fourier  series  applications  are  not  restricted  to 
periodic  functions,  as  explained  in  Gaskell,    special  care 
must  be  taken  in  applying  the  techniques,  the  discussion  of 
which  is  beyond  the  scope  of  this  paper.   Most  texts  on  the 
use  of  the  Fourier  series  provide  the  necessary  instruction 
in  this  area,  e.g.,  Tolstov.  -> 

Although  the  complete  theory  of  Fourier  series  expansions 
are  covered  in  detail  in  such  reference  material  as  Tolstov, 
Gaskell,    Hildebrand,  7  and  others,  a  resume  of  the  simple 
derivation  for  a  periodic  orthogonal  function  is  included 
herein  to  assist  the  reader  in  the  basic  understanding  of  the 
uses  of  the  harmonic  analysis  that  follow.' 

If  we  have  a  function,  such  that  there  is  one  independent 

variable  and  one  dependent  variable,  taken  between  the  limits 

-L  <  0  <  L,  we  could  write  the  function: 

+  L 
z  ■  F(x)     or;     z  ■  S  xdx 

-L 

61 


62 


We  might  illustrate  the  same  function 


A  Z 


-L 
x 


0 


o 


x-0,1,2  .  .  k 


♦L 

x 


k 


The  data  set,  (x.  ,  z  j_)  is  such  that  x   =  -L,  (x.+-,  -  x.  )  is 
constant  and  x^  =  +L,  the  value  of  "k"  being  even. 

We  then  assume  that  an  equation  using  simple  sine  and 
cosine  terms  may  be  written,  which  will  represent  the  func- 
tion, as  follows: 


Z  ■  f(x)  - 


a.- 


I 


nx  2ttx  ntrx 

cos  -V  +  a0  cos  —7 —  +  .  .  .  a  cos  —7 — 
L    2       L  •  n      L 


*  ■>  TJX     x.         2T»X  ^      .    nTTX 

+  b-,  sin  -T-   +  b0  sin  -r —  +  .  .  .  b   sin  — r— 

1.  Lt  c.  Li  n         Li 


or; 


Z  =  f(x)  - 


a 


+  Z 


n-1 


(~  nnx  +  .    .   rrnxx 

(a  cos  —7 —  +  b  sin  -7 — ) 

n      L     n      L 


(1) 


The  degree  to  which  this  equation,  so  written,  will  converge 
upon  the  function,  may  be  shown  to  be  dependent  upon  two 
factors.   The  first  factor  is  the  number  of  harmonics  (n)  to 
which  the  series  is  expanded,  and  the  second  is  the  number  of 
data  points  observed  in  relation  to  the  basic  length.   The 
author  has  achieved  an  accurate  fit  to  the  degree  of  .97^  on 
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an  eighty  .point  function  of  length  2  tt  ,  with  data  values 
ranging  from  -1  to  +1.   The  series  was  expanded  to  forty 
harmonics.   By  holding  either  of  the  stated  factors  constant, 
and  decreasing  'the  other,  a  loss  in  accuracy  of  representation 
was  realized.   The  reason  for  the  "2"  in  the  denominator  of 
the  first  term  of  equation  (1)  will  become  apparent  as  we 
find  the  coefficients  a  and  b   of  the  series.   By  multiplying 
each  side  of  the  above  equation  by  cos  nx ,  and  integrating 
between  the  limits  of  +L  and  -L,  Gaskell  shows  that  the 
equation  may  be  rewritten: 

/   f (x)  cos  ^—-   3x  -  a     f        cos2  ^  dx  -  La 
_L  L        n   -L     '   L         n 


or: 


♦L 


ntrx 


a   -  1  /   f(x)  cos  —^  dx  (2) 

n    L   -L  L 

With  the  a  term  as  the  first,  term  of  the  equation  (l), 

utilizing  the  same  steps,  we  discover: 

♦L  ,     x_        +L  +L  , 

/   f(x)  cos  ^°T;  x  dx  =  /   f(x)  dx  =  /   !&-   dx  ft  La^ 
-L  L         -L  -L  Z   °     ,    ° 

a   -  f  /    f(x)  dx  =  f  /   f(x)  cos  iaiH2L  dx      (3) 
0   L   -L  L   -L  L 

a 
Gaskell  further  notes  that  -—   is  the  average  value  of 

f(x)  over  the  entire  interval.   Using  the  trapezoidal  rule  to 

numerically  integrate  Equation  (2)  over  its  limits,  and 
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2L 


utilizing  equally  spaced  points  (Ax  "  tt~)   where  x  ~  1,  2,  3» 
.  .  .  k,  we  can  approximate  Eauation  (2)  as  follows: 


a   =  =-  •  Ax 
n    L 


o 


'k  cos 


nn 


i=k-l 


+  Y 
1-1 


mix. 

Z.   COS  1 

l      — ; 


z^  +  a,     nn   1vk_1  z 
o k  cos    *  T 

2  i  =1 


nnx. 

i  COS    _J^ 
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(4) 


It  is  shown  by  Whit  taker  and  Robinson,    that  use  of  the 
trapezoidal  rule  produces  coefficients  which  satisfy  the 
least  squares  criterion  of  a  "best"  fit  to  the  existing  data, 
that  is, 


i-k-1 

£    (z.  -  Px. ) 
i     i 
i  -o 


Minimum 


(5) 


where  Px.  is  the  polynomial  approximation  to  be  observed  z.  . 

The  same  reasoning  as  above  can  be  enlisted  to  determine 
the  coefficient  b  ,  except  that  the  sine  terms  are  zero  at 
the  end  points,  (sine  (0)  ■  0)  and  b  becomes: 


,   _    1:k-1      .   nnx. 

b   -  r-  E     z.  sin  l 

n    k  .  _-.     i     — j — 


(6) 


Now,  assuming  we  have  discrete  data  z.  for  each  value 
of  x.  taken  at  equally  spaced  points  between  the  limits,  we 
can  now  compute  the  values  of  the  coefficients,  and  having 
so  done,  substitute  back  into  Equation  (1),  which  becomes  a 
one-equation  representation  of  the  function.   The  coefficients, 
an  and  b  ,  serve  to  characterize,  in  one  way,  the  function. 
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However,  this  means  of  characterization  has  the  distinct 

disadvantage  of  being  dependent  upon  the  phase  angle  of  the 

70 
sinusoids  into  which  the  function  is  resolved.     This 

disadvantage  may  be  overcome  in  several  ways.   Some  of  the 

more  complex  characterization  techniques,  which  we  shall 

discuss  and  compare  in  latter  sections,  are  the  use  of  the 

autocorrelation  function,  and  the  power  spectrum.   These 

functions  are  non-dependent  of  phase  angle,  as  we  will 

later  show.   A  less  complex  technique,  which  still  dismisses 

the  dependency  upon  the  phase  angle  is  discussed  in  a  paper 

on  wind  generated  sea  waves,  distributed  by  New  York 

University.'1    An  estimate,  of  the  power  spectrum  is  made 

in  the  following  manner,  using  the  Fourier  coefficients: 

2      2      2 
Estimated  Power  Spectrum   c    =  a   +  b      (7) 

Although  this  estimate  is  a  characterization  of  the 

function,  and  floes  dismiss  the  dependency  upon  the  phase 

2 
angle,  the  same  reference  notes  that  the  quantity  cn   is  a 

very  unstable  estimate  of  the  energy  at  any  frequency.  ;  Only 

through  a  weighting  process  of  the  c   values,  would  an 

accurate  spectrum,  such  as  that  determined  by  the  Tukey 

method'^  be  recovered.   Lee'-^  uses  a  similar  estimate  of  the 

power  spectrum,  with  the  exception  that  he  divides  the  above 

estimate  by  a  factor  of  four. 

2     2     2 

c„   -  a  ^  +  b  .  /on 

n     n n  ( 8 ) 


^ 


In  that  the  basis  for  our  discussion  dealing  with  auto- 
correlation and  power  spectrum  analysis  evolves  from  this 
same  reference,  and  since  the  division  does  not  affect  the 
term  for  our  purposes,  we  have  used  the  "Lee"  estimate  in 
our  work.   This  usage  will  allow  a  means  of  comparing  this 
simple  estimate  of  the  power  spectrum  with  the  values 
computed  by  more  complex  techniques. 


APPENDIX  B 


DOUBLE  FOURIER  SERIES 


As  in  the  case  of  a  function  with  one  independent 
variable  discussed  in  an  earlier  section,  it  is  also  possi- 
ble to  represent  with  a  Fourier  series,  a  dependent 
variable  as  the  function  of  two  independent  variables. 

z  =  f(x,y)  (1) 

A  rigorous  discussion  regarding  orthogonal  systems  in  two 
variables-  may  be  found  in  Tolstov.^ 

Thus,  a  parameter  such  as  z  may  be  considered  to  be 
oscillatory  in  two  mutually  perpendicular  directions,  such 
as  east-west  and  north-sbuth,  allowing  Equation  (1)  to  • 
represent  the  areal  variation  in  z.   If  we  consider  a 
function  z  to  have  a  fundamental  period  of  2L  along  the  x 
coordinate,  and  2H  along  the  y  axis,  then  following  the  same 
reasoning  discussed  earlier  in  the  single  dimension  case,  we 
can  write  Double  Fourier  series  function: 


m=M    n=N 
P  -      T        E        X 

x'y        m=o  n=o        m'n 


a  cos  nmX   cos  Tiny      . 

IT  H 


b  sin  nmx   cos   Tiny    +    c  cos   Timx    sin   Tiny 


m  ,n 


L 


H 


m  ,n 


L 


Ii 


d  sin   timx    sin  Tiny 

m     n  •*- 


m  ,n 


H 


(2) 
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whe 

re 

X 
m 

,n 

= 

1/4 

X 
m 

,n' 

= 

1/2 

X 
m 

»n 

= 

1     • 
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ra  -  n  _  o 

m  =  o,n>o  or  m>o,n=o 

m>o,n>o 

M,N  number  of  harmonics  considered 

(Should  be  equal  to  or  less  than  one 
half  the  number  of  data  points  con- 
sidered in  the  x  and  y  directions 
respectively.)  ?5 
The  reasoning  behind  the  quarter  and  half  terms  follows  also 
from  the  one  dimension  analysis,  in  that  the  end  terms  con- 
sidered in  a  numerical  integration  are  averaged,  according 
to  the  theory  of  the  trapezoidal  rule.'"   In  the  two  dimen-s 
sion  case,  the  averaging  process  must  be  carried  in  each 
direction,  thus  the  quarter  terms.   Again,  the  problem  is  to 
evaluate  the  coefficients.   Integrating  Eauation  (2),  the 
coefficients  are  determinable  from  the  following  relations: 

Vn  ■  L  s[l  '[[   f<*^>  cos  VJr  cos  T   **«*  (3» 


-n       -L 


cm,n  "  Eh  ;+"  \     «*•*>  oos  ™  sl"  T    W*  (5> 

-n        -Li 


'm.n    "  LH    !Z   '^      «*,*>■  ■"»  f  Sin  =P  d,tt  (6) 


H 
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It  is  possible  to  numerically  integrate  Ecjuations  3-6 
utilizing  the  trapezoidal  rule.   Again  it  must  be  remembered 
that  the  integration  must  be  carried  out  in  two  dimensions, 
and  special  care  must  be  taken  with  regard   to  the  corner 
points.   In  order  to  clari'fy  the  following  derivation, 
Figure  1  displays  the  nomenclature  to  be  used. 


r-* 


J=o 


Z(i,j) 


-4 1— I 1 >    I    4    I 


i-o 


xi 


j=o,l,  2,  ...  .  Jt 
imo,l,  2,  .  .  .  k 


i=k 


Figure  1 


According  to  the  proper  use  of  the  trapezoidal  rule, 
it  can  be  shown  that' for  the  one  dimension  case  (one  . 
independent  variable) 


/    z(x)  dx   ■  Ax 
-L 

a     2L 

where  Ax  -  t— 


^   +  *  i=k-l 

fo Ik  "»   S 

2      1-1 


(7) 


Considering  the  second  dimension  of  depth  in  our  picture, 
when  we  are  concerned  with  but  one  direction,  we  are  able  to 
write  Equation  (7) 


+L 
J*     z(x  )  dx 
-L 


\ 


Ax 


!oJ   zk..i  ■  +  r 

2  1-1 


■i.3 
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(8) 


However,  in  a  specific  case,  where  it  is  required  that 
integration  be  applied  in  two  dimensions,  as  will  be  done  for 
Equation  (3)»  the  following  steps  are  taken: 


1    1    *H   +L 
am,n  =  L  *  H   '  -H  '  _L   zU'y)  dXdy 


»  1 


♦H 


fer  ;-h  Ax. 


^j k_,J.  4 


i-k-1 


1-1 


ijj 


dy 


Ax  Ay 

LH 


3-1-1 


'o,o    k,o  +   o,  .   k,l 
2  2 


-l(A) 


2  tJ    k  J 
2 


2 
(B) 


i-k-1 
1-1 


1  •         Z  .   fi 

1,0     l  ,J? 


(C) 


J'-JUl  l-k-l 

Z       I     z.  . 

3-1    l-i    1,J 


(D) 


(9) 


where:   (A)  term  that  averages  the  corners  in  both 
directions . 

(B)  term  that  averages  opposite  vertical  edges 
( less  corners ) . 

(C)  Term  that  averages  opposite  horizontal 
edges  (less  corners). 

(D)  term  that  sums  all  interior  points. 
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The  number  of  points  in  either  direction  should  be  an 
odd  number,  requiring  "k"  and  "  "  to  be  even  numbers.   Carry- 
ing the  procedure  to  completion  and  reinserting  the  trigno- 
metric  terms,  which  were  omitted,  to  simplify  the  explanation, 
we  can  write  the  coefficient  as  follows: 


a    -l.2g.2H 

m,n   ^Wi   k    1 


'0,0  k,o  +   o,  .   k , Jg. 

2  2     cos  mTT  cos  nn 


3-1-1 
E 

3-1 


'o,.1    k,,1   cos  mTT  cos  nny 
2  — H ' 


(B) 


(A) 


i=k-l 

E- 
i-1 


z. 
1 ,0 


A 


mTTx . 
cos    i 


cos  nn 


(C) 


j-je-i  i-k-i 

T  E 

3*1    i=l 


z. 


1,  J 


mnx.      nny . 
cos  1  cos   J  .1 


L 


J 


(D) 


Similarly,  Equations  (^-6)  become: 


m,n    Ik 


No  (A)  term  in  that  on  the  corner 

— >  o 


points,   sin  mnx1 


L 


(10) 


No  B  term  due  to  the  fact  that 
x  and  x,  are  equal  to  L  and  sine 


o 


k 


term  goes  to  zero. 


+  '7 


V  Zi,o^l,i  Bln   ^ 


cos 


nn 


3  "1-1   i"k-l       .   mux,  „   nny. 
z.  .  sin  1  cos    .1 

1 »  J        T 


(C)  + 
(D)" 


E       E 

3-1     i=l 


H 


(11) 


12 


>     ■»  _ 
'm,n   k 


No  (A)  term  because  sine  terms  go  to 
zero. 


E      o,  .1   k,  .1   cos  miT  sin   J  j 
J-l        2  H 


(B) 


No  (C)  term  because  sine  terms  go  to 
zero. 

'J^-1  i:k'1  mux.   .   miy." 

E  Ij     z,  .  cos  i  sm   J  .1 

3"1'  i=l  •    1,J        L        H 


(D) 


(12) 


m,n   k 


No  (A),  (B),  or  (C)  terms  because  sine 
terms  to  to  zero. 


j=J?-l   i=k-l 
j=l    1-1 


z.  .  sin  i  sin    J  j 

1,J      L         H 


(D) 


(13) 
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APPENDIX  C 

COMPUTER  PROGRAM  TO  DETERMINE  COEFFICIENTS  FOR  SINGLE 
FOURIER  SERIES  (ONE  DIMENSION)  FIT  OF  ORTHOGONAL  FUNCTIONS 

Purpose 
The  purpose  of  this  program  is  to  determine  the  Fourier 
coefficients  a  and  b  ,  when  a  function,  f(x),  is  fitted  by 
a  single  or  one-dimensional  Fourier  series.   The  secondary 
purpose  is  to  calculate  an  estimate  of  the  power  spectrum 
of  the  function,  which  in  turn  is  a  characterization  of  the 
function. 

Language 
Fortran  IV  (IBM  70^0  Computer) 

Symbolic  Dictionary 

Variable     S/A*     I/O**  Description 

X         A        I/O      Independent  variable  data 

point . 

Z  A        I/O      Dependent  variable  data 

point. 

MZ         A         I       Fixed  point  value  for  Z. 

A         A         0       Fourier  coefficient. 

B  A         0       Fourier   coefficient. 


*  S--Single  Variable}   A--Array  of  Variables 
**  I--Input  Variable,;  0--0utput  Variable 
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Variable 


S/A 

A' 


KKK  S 

ZSUM  S 

SUMCS  S 

SUMSIN  S 

SUMCOS  S 


\ 


I/O 
0 


0 


0 


cc 

S 

- 

K 

S 

I 

EL 

s  • 

I 

KOW 

s 

I 

KOO 

s 

I 

PIE 

s 

I 

AZERO 

S 

0 

I 

S 

0 

PMAX 

S 

0 

PMIN 

S 

0 

Description 

Estimate  of  the  function's 
power  spectrum  based  upon  the 
equation : 

2      2      2 
C    -  A   ♦  B 

n     n     n  • 

Estimate  of  input  function 
Z  at  each  data  point,  created 
by  using  calculated 
coefficients. 

Residual  or  difference  between 
input  function  Z  and  estimated* 
function  P,  at  each  data  point. 

A(N)2  •*■  B(N)2 

Number  of  data  points. 

Half-length  of  sample. 

Number  of  samples. 

Sample  index. 

Mathematical  expression  input 
as  being  equal  to  3.^^159265. 

Number  of  harmonics  to  which 
the  Fourier  series  is  to  be 
expanded . 

Sum  of  Z  data  point  values. 

Sum  of  Sine  and  Cosine  terms.' 

Sum  of  Sine  terms  in  equation 
to  determine  coefficients. 

Sum  of  Cosine  terms  in  equation 
to  determine  coefficients. 

First  term  of  Fourier  expansion 

Data  point  index. 

Maximum  value  of  variable  P. 

Minimum  value  of  variable  P. 


Variable 

S/A 

I/O 

ZMAX 

S 

0 

ZMIN 

s 

0 

RMAX 

s 

0 

RMIN 

'S 

0 

N 

s 

am  — 
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Description 
Maximum  value  of  variable  Z. 
Minimum  value  of  variable  Z. 
Maximum  value  of  variable  R. 
Minimum  value  of  variable  R. 
Harmonic  index. 

Program  Routine 
The  program  utilizes  data  points  equally  spaced  along 
the  function  from  a  set  origia,  and  by  using  Fortran  equi- 
valents of  Equations- [k)    and  (6)   from  Appendix  A,  generates 
a  Fourier  series  estimate  of  the  function.   The  Fourier 
coefficients  generated  are  presented  as  characterizing  the 
function,  and  by  combining  them  according  to  the  "Lee" 
estimate,  C('N)2  -  (A(N)2  +  B(N)2)/4,  a  characterizing  power 
spectrum  which  is  independent  of  phase  angle,  is"  determined. 
When  more  than  "one  sample  is  run,  the  C(N)  values  are  grouped 
statistically  for  analysis.   The  program  checks  the  values 
of  the  coefficients  by  creating  an  estimate  of  the  original 
function,  using  Equation   (1)  from  Appendix  A.   A  plot  of  the 
original  function,  the   created  estimate,  and  the  residuals, 
or  differences  between  the  two  is  presented. 
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c 

C  PROGRAM   -FOURIER-          AUGUST  1965 

C  SINGLE  FOURIER  SERIES  PROGRAM  (ONE  DIMENSION) 

C  PROGRAMMED  BY  WD  ALDENDFRFER 

C  ONE  DIMENSION  FIT  OF  3MM  BEADS  DATA___. 

C  400  POINTS     50  HARMONICS 

C  _             


DIMENSION  X(801 ) #2(801 ) *MZ( 801 ) »A(30) »B( 30) ♦C(30»30) tP(801) ♦ 
1IDSH(32 J  .  TITLE(13)»LO(30) »LC(30) .LR(30) »R(801) 
CALL  CLOCK   (0) 


READ  (5,102)  (TITLE( I )»I*1»13) 
102  FORMAT  (13A6)  

READ  (5,1)  K»  EL 
1  FORMAT  (  I10»10X.F10.0)  - 

WRITE(6,6000) 
6000  FORMAT ( 1H1 ♦20X»21HPOWER  SPECTRUM  VALUES) 

KOW  =  16 
<OW  IS  THE  NUMBER  OF  SAMPLES  BEING  TESTED. 

•COO  =  1 

PIE  =  3.14159265 

AK  =  K 

MK  =  K-l  • 

AMK  =  MK 

KKK  =  50 

KKK  IS  THE  NUMBER  OF  HARMONICS 

DELX  =  PIE/( (AK-1. )/2. ) 
555  X(l )  =  -(PIE) 

ZSUM  =  0.0 

SUMCS  =  0.0 

2(1)  =  0.0  

DO  203  III  =  2.K 

Z(  III  )  =  0.0 


X(  I  II  )  =  X( III-l )  +  DELX 

203  CONTINUE 

READ(5,204)  ( MZ ( I ) ♦ I -1 . K ) 

204  FORMAT (80  1 1/801 1/801 1/80  1 1/8011 /ID 
DO  205  JJ  =  1»K 

IF  (MZ(JJ)  -  0  )800»800»801 

800  Z(JJ)  =  -1.0 
GO  TO  205 

801  Z(JJ)  =  MZ(JJ) 

205  CONTINUE 

DO  12  II  -    2.MK 
ZSUM  =  ZSUM  +  Z(  II  ) 
12  CONTINUE 

SUME  =(Z(1)  +  Z(K)  )/2. 

AZERO  =  (2. /AMK)  *(ZSUM  +  SUME) 

DO  4  N  ~    1 »KKK 


SUMS  IN  =  0.0 

SOMF  =  0.0 

SUMCOS  =0.0 

FN  =  N 

SUM=  0.0 

ALGEBRAIC  INDEX  EQUALS  FORTRAN  INDEX  LESS  ONE.  . 

DO  5  I  =  2»MK 

SUMCOS  *  SUMCOS  +    Z ( I )*COS( (FN»  PIE  *  X(I))/EL) 
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SUM5IN  =  SUMSIN  +    Z ( I ) *S IN ( ( FN*  PIE  *  X(I))/EL) 

5  CONTINUE 

SOME  =  ((2(1)  +  Z(K))/2.0)»  COS(FN*PIE) 
c      

C  

C  COMPUTING  THE  VALUE  OF  THE  POWER  SPECTRUM  ESTIMATE  • 

C  ESTIMATE  =  A(N)  SQUARED  +  B(N)_SQUARED  DIVIDED  BY  FOUR 
C 


A(N)  =  2./AMK  *  (SUMCOS  +  SOME) 

B(N)  =  2.0/AMK  *  SUMSIN 

CC  =  A(N)*  A(N)  +  B(N)  *  B(N) 

C(KOO»N)  -  CC/4.0 
4  CONTINUE 

WRITE (6 ,687) (C(KOO,N) ,N=1,50) 
687  FORMAT( (1X,10F10.4) ) 

WRITE(6,6005) 
6005  FORMAT(1HO»11HNEXT  SAMPLE) 

GO  TO  686      

2  NMAX  =  (K-l)/2 

DO  8  I  =  !♦< 


SUMCS  =  AZERO  /  2. 

DO  9  N=l ,KKK 

FN  =  N 

SUMCS  =  SUMCS  +  A(N)  *  COSUFN*  PIE  *  XUM/EU  + 


1B(N)  *  SIN ((FN  *  PIE  *  X(I))/EL) 
9  CONTINUE 

P(  I)  =  'SUMCS 

R(  I  )  =  P(I )  -  Z( I ) 
8  CONTINUE 

GO  TO  213 


C  PLOTTING  ROUTINE 

C 


214  ISTAR  =  -13702990896 . 

IEYE   =   27661634608 
IBLANK  =  -17997958192 
I.DASH   =  -818089008 
DO  10  I  =  1,32 
10  IDSH( I )=  IDASH 

WRITE  (6,103) (TITLE( I ),I=1 tl3> 

103  FORMAT  (  1H1 ,1 3A6// ♦ 17X ♦ 8HOBSERVED »32X , lOHCALCULATED, 
132X,8HRESIDUAL) 

WRITE  (6fl04)  (IDSH( I NI=1»32) ,<IDSH( I  )Yl  =  l,32)7~ 
1( IDSH{ I ) ,1=1,32) 

104  FORMAT ( 5X, 32 A1,9X,32A1,9X, 32 Al) 
ZMAX  =  Z(l) 


ZMIN  =  Z(l) 

PMAX  =  P(l ) 

PMIN  =  P(l ) 

RMAX  =  R(l ) 

RMIN  =  R(l ) 

DO  20  I  =  2»  K  \ 

IF  (Z(I).GT.  ZMAX)  ZMAX  =  Z(I) 

IF  (Z(I).LT.  ZMIN)  ZMIN  =  Z(I) 


-78 


IF  (P(I).GT.  PMAX)  PMAX  =  P(I) 

IF  (P(I).LT.  PMIN)  PMIN  =  P(I) 

IF  (R(I).GT.  RMAX)  RMAX  =  R(I) 
20  IF  (R(I).LT.  RMIN)  RMIN  =  R(I) 

WRITE (6 ,902)   ZMAX »ZMIN »PMAX »PM I N »RMAX »RMIN 

902  FORMAT(1X.4HCHE2»6F10.4) 

DO  250  I  =  lt30 

LO(I)  =  IBLANK 

LC(I)  -  IBLANK 
250  LR(  I  )  =  IBLANK 

KHAF  =  (K  +  1 )/2 

DO  30  1=  lfK 

JZ  =( ( (Z( I )-ZMlN)/(ZMAX-ZMIN) )*  29.)  +  1.   ._ 

JO  =( ( (P(I )-PMIN)/(PMAX-PMIN) )*  29. )  +1. 

JR  =  ( (  (R(I )-RMIN)/(RMAX-RMIN)  )*  29.)  +1. 

LO  (JZ)  =ISTAR 

LC  (JO)  =ISTAR. 

LR  (JR)  =ISTAR 

KK  =  I  -  KHAF 

WRITE(6,105)  KK»(LO(JJ) ♦ JJ=1 ♦ 30 >V  KKVfLC ( JJ ) »JJ*lt30) » 
1KK»(LR( JJ) »JJ=lt30) 
105  FORMAT  (1X»I4>  1HI»  30A1,1HI,  5X  ♦  U  ♦  1HI  »30Al  •  1HI  »  5X  • 
1I4»1HI»30A1»1HI    ) 

LO(JZ)=  IBLANK 

LC(JO)=  IBLANK 

LR(JR)=  IBLANK 
30  CONTINUE 


C      _ u.. 

C  END  OF  PLOTTING  ROUTINE 

213  DO  699  I  =  1»K  

WRITE (6 ,698)   I»Z(I)f  P(I)YrYI> 

698  FORMAT(10X»I5.3F20.8) 

699  CONTINUE 

686  KOO  =  KOO  +  1 

IF(KOO.GT.KOW)   GO.TO  400 

GO  TO  555         _\_  . 

400  CALL  CLOCK (1 ) 

C ._•-—--. .. 

c 

C       GROUPING  OF  C(N)  FACTORS  (POWER  SPECTRUM)  FOR  FREQUENCY 

C       STUDY.  POWER  SPECTRUM  VALUE  =  A*A  +'"B*B  DIVIDED  BY 

C       4.0  AS  DEFINED  IN  THE  TEXT  BY  LEE.  THIS  IS  ONLY  AN 

C       ESTIMATE  OF  THE  TRUE  SPECTRUM. 

C 

C  -----~-z.--~-~^-jr-~~ 


WRITE(6,577) 
577  F0RMAT(1H1,37HFRE0UENCY  DISTRIBUTION  OF  C(N)  VALUES) 

WRITE(6,575) 
575  F0RMAT(lH0t3X»lHN»6X»4HINTl»6X»4HINT2?6X»4HINT3t6X» 


14HINT4»6X»4HINT5»6X»4HINT6»6X»4HINT7»6X»4HINT8»6X»4HINT9» 
26X,5HINT10«5X«5HINT11»5X,5HINT12) 
WRITE(6t576) 
576  FORMAT(1HO»9X»40HLT  .001    .001-. 01   .01-. 02    .02-. 03   "V 
150H.03-.04    .04-. 05    .05-. 06    .06-. 07    .07-.08   •» 
225H08-.09    .09-. 1       GT . 1 ) 
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DO  500  N  *  1,KKK 

TOT  =  0,0 
TOTS  =  0.0 

INT1  =  0 

int?  =  o  : __•_ 

INT3  =  0 

INT4  =  0 . ; ; ■__ 

INT5  =  0 

INT6  =  0  „. 

INT7  =0  •        . 

INT8  =  0  , L 

INT9  =  0 

I  NT  10  =  0  r__: 

INTII  =  0 

INT12  =  0 

DO  501  KOO  =  l»KOW 

COFFF  =  C(KOO.N) 

TOT  ■  TOT  +  COEFF 

TOTS  =  TOTS  +  KOEFF*COEFF)        


IF(COEFF.GT. .00099995)   GO  TO  505 

INT1  =  INT1  +  1 

GO  TO  501 

505  IF(COEFF.GT.. 00995)   GO  TO  506 
INT2  =  INT2  +  1 

GO  TO  501  

506  fF (COEFF .GT ••019995)   GO  TO  507 
INT3  =  INT3  +  1     _:_ 

GO  TO  501 

507  IFKOEFF.GT. .029995)       GO    TO    508 
INT4    =    INT4~+"l" 

_GO    TO    501  

508  IFKOEFF.GT.. 039995)   GO  TO  509 
__j.NT_5  =  IN_T5__+__1 

GO  TO  501 

509  IF  (COEFF. GT..049995J__  GO__TO_510_ 
INT6  =  INT6  +  1 

GO  TO  501 


510  IFKOEFF.GT.. 059995)  GO  TO  511 
INT7  =  INT7  +  1 

GO  TO  501 

511  IFKOEFF.GT..  069995)  GO  TO  512 
INT  8  =  INT8  +  1 

_GO  TO  501  

512  IFKOEFF.GT. .079995)    GO    TO    513 
INT9    =    INT9    +    1 

GO    TO    5,01 

513  IFKOEFF.GT.. 089995)  GO  TO  514 
INT10  =  INT10  +  1 

GO  TO  501 

515 


514 

515 
501 

IFKOEFF.GT.. 1)   GO  TO 

INTII  =  INTII  +  1 

GO  TO  501 

INT12  =  INT12  +  1 

CONTINUE 

AKOW  ■  KOW 

AVE  =  TOT/AKOW 

80 


AVF2  =  TOTS/AKOW 

AVE3  =  AVE*AVE 

SDEV  =  SQRT( ( AVE2-AVE3) ) 


WRITE (6  ,502)  N» I  NT  1 ♦ I  NT  2 ♦ I  NT  3 ♦ INT4» INT5 » INT6 ♦ INT7, INT8» 
1INT9»INT10»INT11»INT12 
502  FORMAT(1HO»I4. 12110) 

WRITE(6,602)      AVE*  5DEV 


602  FORMAT(1HO»7HMEAN  =  »F20.8 * 10X ♦ 16HSTANDARD  DEV.  =  »F20.8) 

c .__. 

C  END  OF  GROUPING  ROUTINE 

c      

500  CONTINUE 

CALL  CLOCK(?)  _ . 

CALL  EXIT 


END 


SENTRY 


\ 


\ 
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SINUS  FOURIER  3ERIE3  PROGRAM 


3TART 


READ  i  TITLE 


T 


READi  K,EL 


KOW  -- 

i  No.    of   3srnpl«g 

■ 

KOO   =  1 

pis  =  5.H159265 


AK  =    K 

MK  -   K-l 

AMK  =    MX 

.     .:. 

. 

KKX  -  No.  of  Harmonics 


OSLX  =  PIE/(K-l)/2 


I *  X(l)  ■  -PIS 


555  fro. 
End  of 

Progrso. 


ZSUM  =  0. 

3UH03  -  0. 

Z(l)-  0. 


205 

III  =  2,X 


Z(III)  ■  0. 


X(III)  =  X(III-l)+DSlX 


«T] 


& 
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j(jj)=  -i.~| 


3UMB  =  (Z(1)»  Z(K))*  } 
AZKRO  =  (2./AMX)  •  (  ZSUH  ♦  3UW) 


N  -  l.KKK 


3UM3IN  »  0. 

3UM003  -  0. 

30MK  -  0. 
FN  «  N 
3UM  ^  0. 


I  =  2,W 


3UM003  -  3UM0O3  +  Z(I)  •  003 
((FT)  •  PIS  •  X(I))/EL) 

3UK3IN  »  3UM3IN  -t-  Z(I)  •  3IN 
((fK  •  PIB  •  X(I))/SL) 


0 
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30ME   -  ((2(1)  +  Z(K))/2  •   U03(FN   ♦   PIS) 


_L 


1 


A(N)  ^    2. /MIX.   *    (3UM003   *•  30ME) 
H(N)  =    2./W1K   •   (3UKSIH) 

I 


A(N)   •   *(K)  -t-B(N)   •   B(N)  (X) 

r  t 


3UM03  -    3UM03  +  A(N)   • 
003((FN  •   PIB  •   X(I))/ 
BL)  -I-  B(N)   ♦  3lK(fN  • 
PIS   •   X(I))/EL) 


\|        H(l)=  P(I)^Z(l)| 


\ 
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3IN')LB  roUAISR  3aRIE.5  i  LOTTIUS  ROUTINE 
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From  M*ln  Fros;rim 


I3TAR   '  1J702990896 
IEYK    •  276616J4608 
IBLANK  =   -1/997958192 
1DA3H  8  -818089008 


/        10 
*\  I  •  1.52/ 

\ 


ID3H(I)  *  IDASH 

I 


Hm 


WHITS  I    HEADING 


WRITEi    I0A3H  LINES 


ZHAX  -  Z(l) 
ZMIN  =  Z(l) 
PHAX  -  P(l) 
PMIN  "  P(l) 
RMAX  :  R(l) 
RMIN   =  R(l) 


Z(I)iZMAX> — ^JZMAX  *  Z(I) 


© 


ZHIN  *  Z(I) 


> 


J 


P(I)iPMAX> — *JP>UX  -   P(I) 

< 


© 
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<^P(I):tnIlT>-         — jifKIK  "  PQ]| 


tmirsi 

z.:ax,  zmin, 

pmax,  pmin, 

rhax,  rmin 


LO(I)  '  I BLANK 
10(1)  -  I  BLANK 
LR(I)    ■  I BLANK 


250 


KHAF^    K  fl/2 


50 
I  -   1,K 


1 


JZ 

-  (((Z(I)-ZHIN)/ 

(ZMAX-ZMIN))    •   29.)+"  »• 

JO 

.   (((P(I)-PMIN)/ 

(PMAX-PMIN))      •    29.)+ 1. 

JR 

-   (((R(I)-RMIN)/ 

(RKAX-RMIN))   •    29.)  +  1. 

"    ' 

LO(JZ)  -  I3TAR 
LO(JO)  -  ISTAR 
LR(JR)  -  I3IAR 


0. 
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KK        I-KHAP 


LO(JZ)  IfLANX 
LO(JO)  IRUNK 
LC(JR)        IHtAKX 


© 


I 


\ 


88 


COMPUTER  PROGRAM  TO  DETERMINE  COEFFICIENTS 
FOR  DOUBLE  FOURIER  SERIES  (TWO  DIMENSION) 
APPROXIMATION  OF  A  SURFACE  AREA 


Purpose 
The  purpose  of  this  program  is  to  approximate  a  surface 
area,  flx,y),  using  a  double  Fourier  series,  and  to  present 
the  Fourier  coefficients  which  are  a  characterization  of  the 
surface. 

Language 
Fortran  IV  (IBM  7040  Computer) 

Symbolic  Dictionary 

Variable     S/A*      I/O**  Description 

X  A         I       Independent  variable  data 

point  along  x-axis. 

Y  A         I       Independent  variable  data 

point  along  y-axis. 

Z  A         I       Dependent  variable  data  point. 

EL         S         I       One-half  of  fundamental  per- 
iod along  the  x-axis. 

H  S         I       One-half  of  fundamental  per- 

iod along  the  y-axis. 

KK         S         I       Number  of  data  points  in 

x-direction. 


* 


S--Single  Variable;  A--Array  of  Variables 
**  I--Input  Variable;  0--Output  Variable 
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Variable 

S/A 

I/O 

LL 

S 

\ 

I 

MMAX 

and 

S 

I 

NMAX 

M 

S 

0 

N 

S 

0 

PIE 

S 

I 

A 

A 

0 

B 

A 

0 

C 

A 

0 

D 

A 

0 

Description 

Number  of  data  points  in 
y-direction . 

Maximum  number  of  harmonics 
to  which  Fourier  series  is  to 
be  expanded.   Should  equal 
one  another. 

Harmonic  index. 

Harmonic  index . 

Mathematical  expression  input 
as  being  equal  to  3^1^159265. 

Fourier  coefficient. 

Fourier  coefficient. 

Fourier  coefficient. 

Fourier  coefficient. 


Program  Routine 
This  program  utilizes  data  points  equally  spaced  upon 
a  grid,  placed  upon  the  surface  to  be  analyzed.   An-approxi-  • 
mation  of  the  surface  is  created  by  expanding  a  double 
Fourier  series  of  the  type  discussed  in  Appendix  B.   The 
Fourier  coefficients  are  presented  as  characterizations  of 
the  surface. 

Input  Generator 
Immediately  following  the  listing  of  the  double  Fourier 
series  program,  is  a  listing  of  a  simple  data  generator  for 
the  program.   This  program  uses  a  random  set  of  coefficients, 
and  Equation  (2)   of  Appendix  B,  to  create,  in  the  form  of 
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punched  output,  data  points  simulating  a  surface  area,  which 
are  the  appropriate  input  for  the  double  Fourier  program. 
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c  

c 

C  PROGRAM   2D  FOUR  I  AUGUST  1965 

C  DOUBLE  FOURIER  SERIES  (TWO  DIMENSION)    PROGRAM 

C  PROGRAMMED  BY  WD  ALDENDERFER   

C 

c  


DIMENSION  X(25)»Y(31)»Z(25»31)»A(5»5) 
DIMENSION  BC5f5) »C(5»5) *D(5,5) 
CALL  CLOCK(O) 

WRITE  (6*600)  • •_ 

600  FORMAT (1 HI »_9H2D  FOURIER  ANALYSIS)" 

WRITE(6»623) 
623  FORMAT(1HO»20'X»26HDISPLAY  OF  CALCULATED  2-D  »  " 
140HCOEFFICIENTS  WHICH  CHARACTERIZE  FUNCTION) 
WRITE(6»654) 
654  FORMAT(1HO»4X»1HM,4X»1HN,7X,6HA(M,N) * 14X »6HB ( MjN ) ♦ 
114X.6HC(M»N) »14X»6HD(M»N) ) 
READ(5»2)   EL»H»KK»LL 
2  FORMAT(2F20.8»  215) 
C       EL  =  FUNDAMENTAL  LENGTH  IN  THE  X  DIRECTION 
C       H   =  FUNDAMENTAL vLENGTH  IN  THE  Y  DIRECTION 
C       KK  =  NO.  OF  DATA  POINTS  IN  X  DIRECTION. 


C       LL  =  NO.  OF  DATA  POINTS  IN  Y  DIRECTION. 

READ(5»1)  (  (X(  I  NYU)  »Z(  I»J)  »  J»l  »LL)  »  I»l  »KK  ) 

1  FORMAT(3F20.8) 
C       X  =  INDEPENDENT  VARTARLE. 
C       Y  =  INDEPENDENT  VARIABLE/ 
C       Z  =  DEPENDENT  VARIABLE.  __ 

MMAX  =  5 

NMAX  =  5 
C      NMAX  AND  MMAX  ARE  THE  NUMBER  OF  HARMONICS  TOWHTCH  THE 
C      SERIES  IS  TO  BE  EXPANDED. 

PIE  =  3.14159265 

AKK  =  KK-1      

ALL  =  LL-1 

M<K  =  AKK 

MLL  =  ALL 
AMK  =  AKK»ALL 
DO  22  M  =  1»MMAX 

FM  =  M-l  

DO  22  N  =  1»NMAX 
FN  =  N-l 
ABLE  =  0. 
ARAKF  =  0. 
ACHAR  =  0.0 

ADELT  =  0.0 

BCHAR  =  0.0 
BDFLT  =  0.0 
CBAKE  =  0.0 
CDELT  =  0.0 
DDELT  =  0.0 

PIM  a  pie  *  FM  

PIN  =  PIE*FN 
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c 

C  CALCULATE    A(0»0) 

C 


ABLE    =    ((((Z(l«l)    +    Z(KK»l))/2.0)    +    UZ(1*LL)    +    Z ( KK »LL ) ) /2 .0 )  ) 

1/2.0)    *    COS(PIM)    *    COS(PINJ 

DO    20    I     =    2>MKK 

DO    ?0    J    =    2»MLL  _ 

ARGX  =  PIM  *  (X( I )/EL) 
ARGY  =  PIN  *  (Y(J)/H) 
CX  =  COS (ARGX) 


CY    =    COS(ARGY) 

c 

c 
c 

SX    =    SIN(ARGX) 
SY    =    SIN(ARGY) 
IFfI.GT.2)       GO 

CALCULATE    COEF 

TO  400 


BY  PARTS. 


ABAKE  =  ABAKE  +((Z(1*J)  +  Z(KK»J))/2.0)  *CX  *CY 

CBAKE  =  CBAKE  +  ((Z(1»J)  +  Z  (KK » J  )  )  /2. 0  )  *  CX  ♦  SY 

400  IFU.GT.2)   GO  TO  401 

ACHAR  =  ACHAR+  UZ(I»1)  +  Z ( I »LL )  ) /2 .0 )  *CX*CY 
BCHAR  =  BCHAR  +  ((Z(Itl)  +  Z ( I »LL ) )/2 .0 )*SX*CY~ 

401  AOELT  =  ADELT  +  Z(I»J)  *CX*CY 
BDELT  =  BDELT  +  Z (  I  .  J ) *SX*CY 
CDELT  =  CDELT  +  Z(I»J)  *  CX  *  SY 
DDELT  =  DDELT  +  Z(I»J)  *  SX  *  SY 

20  CONTINUE 

A(M#N)  =  (ABLE+  ABAKE+  ACHAR  +  ADELT)  *  "(4*0?  { AKK*ALL )  F 
B(M»N)  =  (BCHAR  +  BDELT)  *(4.0/AMK) 


C(M»N)  *  (CBAKE  +  CDELT)  *  (4.0/AMK) 
D(M*N)  =  (DDELT)  ♦  (4.0/AMK) 
22  CONTINUE 

C___  _  

C      WRITE  VALUES  OF  COEFFICIENTS. 
C 


DO  25  M  «  ltMMAX 

MM  =  M-l 

DO  25  N  =  1»NMAX 

NN  =  N-l  •  __ 

W  R  I T  E  (  6  ♦  2  7  )   MMtNN»A(M»N)tB(M»N)»C(M»N)»D(M»N) 

27  FORMAT(1HO»2I5,4F20.4)  

25  CONTINUE 

CALL  CLOCK(l) 

CALL  EXIT 


END 


SENTRY 
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c 

C  PROGRAM  2D  GENERATOR 

C  PROGRAM  TO  GENERATE  DEPENDENT  VARIABLE  Z  AS  A  FUNCTION 

C  OF  TWO  INDEPENDENT  VARIABLES  X  AND  Y.   PROGRAM 

C  USES  DOUBLE  FOURIER  EQUATION  AND  GIVEN  COEFFICIENTS. 

C  COEFFICIENTS  USED  IN  THIS  PARTICULAR  PROGRAM  ARE 

C  TAKEN  FROM  PRESTON-HARBAUGH  PAPER. 

C  PROGRAMMED    SEPTEMBFR  1965. 

C  • 

DIMENSION  A(5»5)»B(5»5)«C(5»5)»D(5t5)»T(25»31) 

CALL  CLOCK  (0) 


DO  100  N«l»5 

DO  100  M»l.»5 

FN=N 

FM=M 

A(M»N)  =  (.3  +  .02*FN  +  .11*FM)*100. 

B(M»N)  =  (.2  +  .04*FN  +  .08*FM)*100. 

C(M»N)  =  (.1  +  .06*FN  +  .05*FM)*100. 

D(M»N)  =  (.08  +  .05*FN  +  .02*FM)  *100. 

100  CONTINUE 

DO  101  M  alt 5 
C(M»1)  =  0.0 

B(1»M)  ■  0.0 _  

D(1»M)  =  0.0 

101  D(M,1)  =  0.0 


PI  =  3.14159265 

SUMZ  =  0.0 

AH  =  15. 

AL  =  12.  

WRITE(6»984) 
984  F0RMAT(1H1»29HINPUT  VALUES  OF  COEFFICI ENTS. ) _ 

WRITE(6,982) 

982  FORMAT(1HO»5X»81HM     N  A(Y»N)_    B(M»N) 

1       C(M»N)  D(MtN)) 

DO  300  N  =  1»5  ■_ 

DO  300  M  «  1»5 

WRITE(6,981)M»N»A(M»N)tB(M»N)»C(M»N)»D(MtN[ 
981  FORMAT(2X»2I5»4F20.8) 
300  CONTINUE 


PIL=PI/AL 
PIH=PI/AH 
XC  =  -12. 
YC  =  -15. 

XINC  =,1.0 
YINC  =  1.0 


WRITE(6»983) 

983  F0RMAT(1H1»11X»24HX(  I  )  Y(J)  » 

122H  Zff«~J)*l3H        !     J) 

DO  410  I  -    1*25 
DO  405  J  «  It  31 
DO  200  N=l»5 
X  =  XC 

Y  =  YC  

FNT=N-1 

PINH  ■  PIH*FNT 


94 


DO  200'  M  =  l»5 

FMT=M-1 

PIML  -  PIL*FMT 

IF(N-1)140. 120*140 
120  IF(M-1)150»130»150_ 
130  AMDA=.25 

GO  TO  170 
140  IF(M-1 ) 160*150tl60 
150  AMDA=.50 

GO  TO  170 
160  AMDA=1.00 


170  ARGX  =  PIML*X 

ARGY  -    PINH*Y 

CX  a  COS (ARGX) 

CY  =  COS (ARGY) 

SX  =  SIN(ARGX) 
_SY  =  SIN (ARGY)     

ACC  =  A(M»N)»CX*CY 

BSC  =  B(M»N)»SX*CY 

CCS  =  C(M»N)*CX*SY 

DSS  =  D(M.N)*SX*SY 

TS  =  AMDA*(ACC  +  BSC  +  CCS  +  DSS") 

ZC  =  ZC  +  TS 

200  CONTINUE 

WRITE(6,980)  X,Y»ZC»I»J 
980  F0RMAT(1X«3F20.8#2I5) 

WRITE(7,950)  X,Y»ZC_ 
950  FORMAT(3F20.8) 


J_ 


ZC  =  0.0 
405  YC  =  YC  +  YINC  ._ 

YC  =  -15. 
410  XC  =  XC  +  XINC  

T(l»l)  ■(f(ltl)  +  T(l»31)  +  T(25»l)  +  tf25»31))/4.6 

T(l,31  )  =  0.0 


T(25»l)  =  0.0 
T(25»31 )  =  0.0 
DO  700  J  =  2*30 

T(1»J)  =  (T(lf  J)  _+__TJ25»J))/2«0 
T(25»J)  *  0.0 
700  CONTINUE 


DO  701  I  =  2»24 
J(M)  =  (TIM)  +  T(I*31))/2.0 
T(  1,31)  =  0.0 

701  CONTINUE 

DO  702  I  =  1»25 
DO  702  J  =  1»31 


SUMZ  =  SUMZ  +  T( I, J) 
702  CONTINUE 

AOO  =  (4.0*SUMZ )/720. 
WRITE(6»9999)   SUMZ»  AOO 
9999  FORMATdHl  »18HSUM  OF  Z  VALUES  «  »F20.8  »5X  ♦9HA707oT  *"  VF20.8  ) 

CALL  CLOCK  (!) 

490  CALL  EXIT 
END 
SENTRY 


DOUhLS    TO   RISR   jKrlIK3    FH03RAM 
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5TART 

3  . 

|    VRITBi    TITLE 
^    «n<)   HKAOIKOS 

Trsadi  el,  h,  kk,  uT] 

. —      ~r~ —  . 

<    RdADi    X(I),   V(J),    Z(I,J) 

"  1 


NMAX   »  NO.    HARMOt  I'.  J 
MMAX  =   NO.    HARMON  103 


P»-»*| 


ABLE 

■  0. 

ADELT 

*  0. 

A  BAKE 

-  0. 

ACHAR 

-  0. 

BOHAR 

1 

•  0. 

BDELT 

-   0. 

1      OBAKE 

T     0. 

ODSLT 

*  0. 

•    DDELT 

=  0. 

d) 
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PIM  -  PIE  ♦  m 
PIN  =•  PIS  •  FN 


ABLE  -  (({(2(1.1)  *  ?(KK,l))/2.0) 

((Z(1,LL)  +  Z(KK,LL))/2.0))/ 
2.00  *  003  (PIM)  *  003  (PIN) 


20 
I  -  2.KKK  X- 


20     \ 
J  =  2.HLL  )*■ 


AR1X  -  PIK  •  X(I)/BL 
ARGY  =  PIN  ♦  X(I)/H 


X 


OX  ^  C03(  AROX) 

OY  ■  003(  AR3Y) 

3X  »  3IK(  AHOX) 

9Y  =  3IN(  A ROT) 


ABAKE        APAKS      (('(I, J) 
Z(KK,J))/2.0)   •  OX  •  OT 

OBAKS        CHAKB      ((Z(I,J) 
Z(KK,J))/2.0)   •   OX   •  AY 


© 


AOHAR        AOHAR     ((2(1,1) 
t(I,LL))/2.0)  •  OX  •  OY 

BOHAR        MHAR     ((Z(t,l) 
Z(I,LL))/2.0)  •  8X  •  OT 


© 
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© 


Q 


ADSLT  =  ADELT  +7.(I,J)  •  ;«  •  Ul 
-OELT  =  BOSLH-Z(I.J)  ♦  3X  •  OY 
0015LT  -  ODBLT+Z(I,J)  •  OX  •  3Y 
DOKLT  =  D05LT»-Z(1,J)    •     It   «   9Y 


©= 


J 


A(M,N)   =  (AHL3  +  ABAKE-f  AOHAR-t  AD5LT) 
•    (4./AKK    •    ALL) 

H(H,N)  =   (HOHAm-MOBLT)   ♦    {'l./HUK) 

O(M.N)  =   (OHAKS  +  CWiLT)    •    (4./AMX) 

D(H,H)  ■    (OnSLT)   •    (b./AMK) 

0 


0 


._  J 


25 
M  *  l.MMAX 


MM  =  M-l 


25 
N  >  l.NMAX 


NN  =  N-l 


WRITEi  MM,  NN,  A(M,N), 
B(H,K),  0(M,N),  D(M,H) 


I 


OALL  EXIT 
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COMPUTER  PROGRAM  TO  DETERMINE  COVARIANCE 
FUNCTION  AND  POWER  SPECTRUM 


Purpose 
The  purpose  of  this  program  is  to  determine  the 
covariance  function  and  the  power  spectrum  of  a  one- 
dimensional  function,  f(x).   The  secondary  purpose  is  to 
plot  the  power  spectrum  and  the  confidence  bands  of  the 
spectrum,  to  facilitate  analysis  of  the  function. 

■  » 

Language 

Fortran  IV  (IBM  7040  Computer) 

Symbolic  Dictionary 

Variable     S/A*      I/O**  Description 

X  .        A         I       Value  of  the  function,  at  each 

data  point.   Read  in  as  fixed 
point,  mx . 


\ 

Series  index. 

Lag  Index. 

Covariance  estimate. 

Raw  estimate  of  the  power 
spectrum.  - 

U         A         0       Smoothed  estimate  of  the 

power  spectrum. 


I 

S 

- 

K 

.      S 

- 

C 

A 

0 

L 

A 

0 

*  S--Single  Variable;  A--Array  of  Variables 
**  I--Input  Variable:  0--Output  Variables 


\ 


Variable 

S/A 

I/O 

N 

S'  . 

I 

M 

s 

I 

EDF 

A 

I 

TL95 

A 

I 

TL90 


LAG 


UMAX 

UMIN 

SPANA 

SPANB 


A 


TU95 

A 

TU90 

A 

KOW 

S 

KOUNT 

S 

LA 

A 

s 

s 


0 
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Description 

Number  of  data  points. 

Maximum  lag.   (Less  than  N/3) 

Equivalent  degrees  of  freedom. 

Table  of  confidence  limits, 
lower  95  per  cent.  » 

Table  of  confidence  limits, 
lower  90  per  cent. 

Table  of  confidence  limits, 
upper  95  per  cent. 

Table  of  confidence  limits, 
upper  90  per  cent. 

Number  of  samples . 

Sample  index. 

Storage  for  data  identifi- 
cation . 

Degree  of  offset  by  which 
the  function  is  correlated 
with  itself,  to  determine 
covariance  estimate. 

Maximum  value  of  smoothed 
power  spectrum  in  each  fre- 
quency  band. 

Minimum  value  of  power 
spectrum . 

90  per  cent  confidence  band 
of  smoothed  power  spectrum. 

95  per  cent  confidence  band. 


Program  Routine 
This  program  utilizes  data  points  which  represent  the 
value  of  the  function  to  be  analyzed  at  equally  spaced 
intervals  along  the  function.   Using  a  numerical  integration 


100 

techniaue,  the  trapezoidal  rule,  the  values  of  the 
covariance  estimate  are  determined  according  to  Equation  (3«25) 
of  this  paper.   A  rough  estimate  of  the  power  spectrum  is 
computed  using  Equation ( 3.26) and  smoothed  according  to 

Equation  (3*27).  The  table  of  confidence  limits  presented  by 

77 

Granger    is  read  into  the' program  and  the  ninety  and 

ninety-five  per  cent  confidence  bands  are  determined. 

The  smoothed  power  spectrum  with  these  confidence  bands  is 

displayed  graphically. 
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c 

C  PROGRAM   COVAR    FOR  COMPUTING  COVARIANCE  FUNCTION 

C  AND  POWER  SPECTRUM.   PROGRAMMED  SEPT.  1965  BY 

C  FW  PRESTON  AND  WD  ALDENDERFER. 
C 

C  r... 

C  NOMENCLATURE 

C  K  -    LAG  TNDFX 

C  I      =  SERIES  INDFX 

C  J        WEIGHTING  FACTOR  INDEX 

C  C(K)   =  COVARIANCE  ESTIMATE 

C  L(J)   =  RAW  ESTIMATE  OF  POWER  SPECTRUM 

C  U(J)   =  SMOOTHED  ESTIMATE  OF  POWER  SPECTRUM 

C  M      =  MAXIMUM  LAG  (LESS  THAN  N/3) 

C  N        NUMBER  OF  DATA  POINTS 

C  LA(I)  =  STORAGE  FOR  DATA  IDENTIFICATION 

C  EDF(I)  =  EQUIVALENT  DEGREES  OF  FREEDOM 

C  TL95(I)=  TABLE  OF  CONFIDENC  L  IM  I  TS  _^_  LOWER-  95PERCENT_ 

C  TL90(I)=  -  LOWER-  90PERCENT 

C  TU90U)=  _-  UPPER  90  PERCENT 

C  TU95(  I  )=  "-"UPPER"  95  PERCENT 

C  KOW  IS  THE  NUMBER  OF  SAMPLES 

C 

C 
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400 

99 

100 

DIMENSION   X ( 3200 ). C ( 301 ) »L( 302) »U( 301) »EDF ( 13 ) »TL95 ( 13  )  » 
1TU90(13) fTU95(13) »LA(13) ,J0W1(53) , JOW2 ( 53 ) »TL90( 13 ) »MX(3200) 

REAL  L 

CALL  CLOCK   (0) 

KOW  =  2 

KOUNT  =  1  _ 

READ    (5  ♦110)     (EDF(I)»TL95( f WTL90(I )»TU90< I) »TU95fI) »I  =  1»13) 

FORMAT  (  F5.  1  ,4F5..2)_ 

READ  (5,99) (LA( I )Vl  =  l ♦ 13 )  ' 

FORMAT (13A6) 

READ  (5,100)  N»M"" 

FORMAT  (215) 

READ  (5,101)  (MX( I )»I=1»N) 
101  FORMAT(80I1 )  ' 

DO  700  11=  1,N 

X( II)  =  MX( II  )      , 
700  CONTINUE 

DO  10  I  =  1,N 


IF  (X(I)  -  0.  )lltll»lO 
11  X(I)  =  -1.0 
10  CONTINUE 

C  ••_ 

C      COMPUTING  COVARIANCE  ESTIMATE   C(K) 
C 


FM  =  M 

MPl  =  M+l 

FN  =N 

J  =  0 

AK  =  2.0*FN/FM 

DO  910  I  *  It  13 

J  =  J+l 

IF(AK.LT.EDF(I ))   60  TO  320 

910  CONTINUE 
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320  J  =  J-l 

CHG  =  (AK-EDF(  J)  )/(EDF(  J+l  }-EDF(JM 
T90L  =  TL90(J)  +  { TL90 ( J  +  l ) -TL90 ( J ) )*CHG 
T95L  =  TL95(J)  +  ( TL95 ( J+l ) -TL95 ( J ) ) *CHG 
T95U  =  TU95U)  +  ( TU95 ( J+l )-TU95 ( J ) )*CHG 
T90U  =  TU90(J)  +  (TU90( J+l )-TU90( J) )*CHG 

DO  210  KK  -    1.MP1 

<  =  KK-1 

51  =  0,0 
NMK  =  N-K 
FNMK  a  NMK 
KP1   =  K+l 

DO  205  1=1  ♦NMK 

IPK  =  I+K 

205  SI  =  SI  +  X(  I  )*X(  IPK') 

52  =  0,0 

DO  206  I=KP1»N 

206  S2  =  S2  +  X( I  ) 

53  =  0.0 


DO  207  I=1»NMK 
207  S3  =  S3  +  X(I  ) 
210  C(KK)  =  (SI  -  (S2*S3)/FNMK)/FNMK 

C 

C      COMPUTATION  OF  RAW  ESTIMATE  OF,  POWER  SPECTRUM- L  {  J  ) 

C 


MN1=  M-l   .  • 

DO  230  JJ  =  ltMPl 

J=  JJ-1 

FJ=J 

SI  =0.0 

DO  220  K=1»MN1 

F<  =K 
220  SI  =  SI  +  C(K+l)*COS(3.1415926*FK*FJ/FM[ 
230  UJJ+1)  =(2.0*S1  +  C  ( 1  )  +  C(M+l)*COS(3Vl4T5926*FJ)")/6.2831952 

L(1J=L(3) 

L(M+3)  =  L(M+1) 


C .  _____ 

C      COMPUTATION  OF  SMOOTHED  POWER  SPECTRUM  ESTIMATE 
C 

MP2  =  M+  2 

DO    240    J    =2»MP2 
240    U(J)    =    .25*L(J-1)    +  ~.50*L(J>    +~.2l5*CU+lT  ' 

C_ 

C      WRITE  HEADINGS 

c ____         1: 

WRITE(6,103) 

103  FORMATdHl  ,36X»25HP0WER  SPECTRUM  ESTIMATION^ 
137X»24HFOR  DISCRETE  TIME  SERIES/31X* 
237HUSING  TUKEY  HANNING  WEIGHTING  FACTORS///) 

WRITE(6,104)  (LA( I  )  »I  =  1»13) 

104  F0RMAT(1X»13A6///) 
WRITE(6»106) 

C 

C      WRITE  POWER  SPECTRUM  AND  TOLERANCE  ESTIMATES 

C 
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DO  250  I  =2»MP2 

LAG  =  1-2 

U90L  =  T90L*U( I ) 

U95L  ■  T95L*UU) 

U90U  =  T90U»U( I ) 

U95U  =  T95U*U(I) 
250  WRITE(6,105)  LAG»C(  1 1  )  ,L  (  I  )  ♦  U90L^U  U  )  ♦U95U  »U95L  »U  (  I  )  »U90U 
105  FORMAT(1X,I3»8(2X,F10.4) ) 
C 

C  PLOTTING  ROUTINE 

C 

764  ISTAR   =-13702990896 

IMINUS  =  29809118256 

IPLUS   =  17997958192 

IBLANK  =  -17997958192 

II      =  27661634608 

WRITE  (6tl03) 

WRITE  (6,104) (LA( I ) ,1=1 ,13) 

WRITE  (6tl07) 


107  FORMAT (19X»23HSM00THED  POWER  SPECTRUM , 34X i 
123HSM00THED  POWER  SPECTRUM  /15X, 
230HWITH  90  PCT.  CONFIDENCE  LIMITS,27X, 
330HWITH  95  PCT.  CONFIDENCE  LIMITS/6X, 


451H0     1     2     3     4.    5     6     7    ™W~       ~9~        0, 
56X,51H0     12     3     4     5     6     7     8     9     0) 
106  FORMATf 1H0,39X, 

150HSMOOTHED  SPECTRUM  WITH  INDICATED  CONFIDENCE  LIMITS/ 

234X.7H90  PCT. » 17X , 7H90  PCT.,5X,7H95  PCT. ♦ 17X ,7H95  PCT./8X, 

310HCOVARIANCE»5X,3HRAW,9X,5HLOWER»19X,5HUPPER»7X, 

45HL0WER,19X,5HUPPER/1X,3HLAG,5X,8HSPECTRUMV4X,8HSPECTRUM,~ 

56X,5HLIMIT,6X,7HAVERAGE,6X,5HLIMIT,7X,5HLIMIT,6X,7HAVERAGE, 

66X,5HLIMIT) 

C 

C      PLOTTING  OF  SPECTRUM" 

C 


DO  280  1=1,53 

JOW1 ( I )= IMINUS  

280  JOW2( I )=IMINUS 

WRITE(  6,108)     (JOWKI  )  ♦  1  =  1 ,53 ) , ( JOW2 < I) • 1*1 »53 ) 
108    F0RMAT(5X»53A1,5X,53A1 ) 

DO    290    1=2,52 

JOWl (I )=IBLANK 
290  J0W2(I  )=IBLANK   __ 

JOWl (1 )=II 

J0W1(53)=II 

JOW2(l)  =11 

J0W2(53)=II 

UMAX  =  U(2 ) 

UMIN  =  U(2) 


DO  300  I  =  3,MP2 
IF(U( I ).GT.UMAX)  UMAX=U(I) 
IF(U( I J.LT.UMIN)  UMIN=U(I) 
300  CONTINUE 

SPANA  =  AL0G10(UMAX*T90U)  -  ALOGl  0  (  UMI  N*T90'L ) 
SPANB  =  AL0G10(UMAX*T95U)  -  ALOGl 0( UMIN*T95L ) 
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ZA    =    ALOG10(UMIN*T90D  

ZB  =  ALOG10(UMIN*T95L) 

DO  310  I  =2»MP? 

IA1  =  ( (ALOG10(U( I )*T90L)  -  ZA  ) /SPANA ) *50.  +  2.0 

IA2=( (ALOG10(U( I ) )-ZA)/SPANA)*50.  +  2.0 

IA3*( (ALOG10(U( I )*T90U)-  ZA ) /SPANA ) *50 .  +  2.0" 
JOWKIA1)  =  IMINUS 

J0WKIA2)  =ISTAR 

J0WKIA3)  =IPLUS 

IB1=( (ALOG10(U( I )*T95D  -  ZR ) /SPANB )*50.+  2 .0 

IB2  =  ( (ALOG10(U( I ) )-ZB)/SPANA)*50.  +  2.0 

IB3=( (ALOG10(U( I)*T95U)-ZB)/SPANB)*50.  +  2.0 

J0W2(IR1)  =  IMINUS 

JOW2UB2)  =  ISTAR 

JOW2CIB3)  =IPLUS 

LAG  =1-2 

WRITE (6, 109)  LAG»(J0W1( J) ♦ J=l »53 ) »L AG, (J0W2(J) »J=1»53) 
109  FORMAT (1X» 13 »lXY53Al»IX» 13 »1X»53A1) 

J0W1UA1)  =  I  BLANK 

J0WKIA2)  =  I&LANK 

J0WKIA3)  *  IBLANK 

J0W2UB1)  =  IBLANK 

JOW2UB2)  =  IBLANK 

JOW2(IB3)  *  IBLANK 
310  CONTINUE 


C 

c 
c 


167 


765 


END 

UNT.EG 
CLOCK 
=  KOU 

OF  PLOTTING 

ROUTINE 

IF(KO 

f.KOW) 

GO 

TO 

76  5 

CALL 
KOUNT 

(2) 
INT  +  1 

•             t 

GO  TO 
CALL 

400 
EXIT 

END 


PROOJUM  TO  W',TC'RMIHB  OOVARIANOI 
P'JMOTIOK   AND  FOYER   SPEOTWM 
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START 


KOV  •  NO.    SAMFLB9 


KOUNT  «   1 


RKADl   Ei>F(I),   TL95(I),   TL90(I), 
TIJ90(I),  TU95(I) 


© 


1    400  RKADl    TITLE 


READi    N,M 


READi   MX(I) 


700 
II  =  J.H 


X(II)  •   MX(II) 


x(i) 

■  -1. 

-_.J 


~_J 


0 
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j  -  j-i 


OHO  •    (AK  -  EDF(J))/ 

(EDF  (J   10  -  EDF  (J)) 


T90L  =■  TL90(J) 
(TL90(J+1)  -  TL90(J))   •  CHO 

T95L  »   TL95U) 
(TL95(J+1)  -  TL95(J))   *  OHO 

T95U  =  TU95(J) 
(TU9?(J+1)  -  TU95(J))    ♦  OHO 

T90U  =  TU90(J) 
(TU90(J+1)  -  TU90(J))   •  OHO 


\ 


210 
M  =  1.HP1 


1 

' 

1  KK  «  KK-1 

< 

' 

|     31  =  0. 

■ 

i 

HMX  ■  N-I 

— 

1 


I 


\ 


FNHK  »   KMK 

KP1  -    K  *  1 

205 

I  =   l.NHK 


r 


ipk  =  i-rn 


I 


D 


31 


!ll  +  X(I)   *   X(IFK) 

205    I 


/  3o6  V 

\i  -  kpi,  *  r 


32  =  S2  +  X(I) 


207 

i  -  i,roo    *  - 


3J=  3J+  X(I) 


L-i? 


0(KK)  «  (31  -   (32»3J)  /   rHMX  /  FNHI) 


^ 


I 


nm»  m  -1 


230 

JJ  -  1.MP1 


5"  6 


10? 
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JJ  =  J-l  1 

1 

1  w  -  7| 

81  ■    0.0 

' 

/          220        \ 

•  - 

re  -  k 

< 

31   «   31  4-  0(K  -*-  1) 

*  ooiiv  •  n  *  ti/rvL) 

- — 

(  220  V-    —    - 

I . 


L(JJ  ♦  1)   ■    2.0  *   81    *  0(1) 
+  0(M  *  1)    ♦  003(TT»   W)/8  •  TT 


0. 


/ 

UX)  --  L(J) 

Km  ♦■  5)  -  L(M  +•  l) 

' 

HP2  »  m  +•  a 

' 

260 

J  *  2,MP2 


<? 


u(-0  -  .25  *  KJ-1) 

+  .5  •  L(J)  +  .2J  *  K-f  «•  1) 


2*0    >-    —    — 
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U90L  =  T90L  •  U(  I ) 
U95L  =  T9tX  ♦  U(I) 
U90U  =  T90U  ♦  U( I ) 
U95U  =  T95U  •  U(I) 

'          , 

WRITE i  LA 

L(I),  U90L, 

U95L,  U( 

1.  C(I-l) 
U(I),  uwu, 
I),  U95U 

•   / 


(    25°    ) i 1 


PLOTTING  ROUTINE 


IOUHT  ■ 

=  KOUNT 

+  1 

" 

C   00 

TO  i00 

) 

"OOVAR"    PLOTTItO   BOUTINS 


•    110 


\ 


-lURf 


I3TAR  =  -1  '^2990896 
IHINU3  -  29009118256 
IPLU3  =  17^7958192 

I  PLANK  -    -17997953192 

II  =       27'6l6}'l601 


JOWI(I)  --    IMINU3 
J0W2(I)  -   IHINU3 


I 


JOtfl(I)  -    I BLANK 
J0»2(I)  =    IhLANK 


0- 


I 


J0W1(1)  -=  II 

J<wi(5J)=  II 

J0W2(1)=  II 

J<W2(5J)=  II 


UMAX  "  U(2) 

imiN  -   U(2) 


© 


Ill 


I 


\ 


U.IAX  *  U(I) 


T 


(I)i    J.,I^>.   *    >  U1IN  •=  'J(I) 


3PANA  ■=■  ALOOlOi'JMAX   «   T90U) 
-AL0G1C)(UMIN  ♦   T90L) 

3FANB  =  ALOj10(UMAX  ♦  T95U) 
-AL0310((JMIN  •  T95L) 


ZA  «   AL0010(Ut:iN   •  T90L) 

zb  =  alooio(um:ik  •  T95L) 


510 
I  =  2.MP2) — 


"1 

4 


IA1  »  ((ALOG10(U(I)    *   T90L) 
-ZA)/:JPANA)    •    50.  +■  2. 

142  =  ((AL0310(U(I))-ZA 
/3PANA   *   50.  +    2. 

1*5  -  ((AL0010(U(I)    *   T90U) 
-ZA)/3PAKA)    «    50.+  2. 


JOWl(IAl)  -    IMIKU3 
J0fl(lA2)  =   ISTAR 
JOWl(IAJ)  =•  IPLU3 


IB1  =   ((»L0Q10(U(I)    •  T95L) 
-!B)/3PAKB)    •   50.      2. 

IB2  -  ((ALOOIO(0(I))-ZB) 
/JrAKR)   •  50.        2. 

IB>=    ((ALO0l0(U(I)   *  T95U) 
-VK)/SPABB)   •   «0.      2. 


JOO(lhl)  ,    IMINU3 

•'■  *?.[IK>)  r.    I3TAJ» 

JOI«(TtH)  =     1K.US 

I 
® 
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I  WRI  TF.i    l\    ,JC-ll(,i) 
LA1,J0W2(,T) 


^~— 

" 

JOVl(IAl) 

=    I-U'K 

JOWl(IA2) 

=    IhLAKK 

JOtfl(IAJ) 

=   IHLAKK 

v 

r 

J0V2(IH1) 

=    It-LAMt 

J0W?(IR?) 

-    I'-LAfK 

JOU?(IB?) 

S   InLAKK 

<? 


_  J 
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APPENDIX   D 

PICTURES 


iOUki 


Sample. X  Three  Millimeter  Beads 
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Saraple  Y  Two  Size  (4mm.  and  5mm.)  Beads 
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Sample  Z  Three  Size  (3mm.,  4mm.,  and  5mm.)  Beads 
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1""         '  "w  i^tvm-  t\\vwi.mvw  ^y.^i.T-ym.^.^y.  ^^  ■&-.*<$ 


V 
v.,     ' 


■ 


^.v<#~"~»% 


7*y>«w:*& 


J        *■-.   '-I'-  ■{■  ■    ■  ■      ■-      ■■■ —        |f.->,lf,  :■„,  ■,>fc!ll,-^-         -  jL~^.J 


Photomicrograph  of  porous  medium 
(sandstone)  92X.  Courtesy  of  Humble  Oil  and  Re-* 
fining  Company.  Porosity  26.3  per  cent.  The  line 
is  0.17  cm  in  length. 


^ 
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POWER  SI'ECIRUM  ESTIMATION 

FOR  DISCRETE  TIME  SERIES 

USING  TUKEY-HANNING  WEIGHTING  FACTORS 


TWO  SIZE  BEADS 

a) 


40  ROWS  MIXED 


LAG 


40 


SMOOTHED  SPECTRl 

)M  WITH  INOICATEO  CONF 

IDENCE  LIMITS 

90  PCT. 

90  PCT. 

95  PCI. 

95  PCT. 

COVARIANCE 

RAW 

LOWER 

UPPER 

LOWER 

UPPER 

LAG 

SPECTRUM 

SPECTRUM 

LIMIT 

AVERAGE 

LIMIT 

LIMIT 

AVERAGE 

LIMIT 

0 

0.9409 

0.7021 

0.6647 

0.7353 

0.79U 

0.6367 

0.7353 

0.7970 

1 

0.6656 

0.7684 

0.6291 

0.6959 

0.7488 

0.6027 

0.6959 

0.7544 

2 

0.4372 

0.5449 

0.5591 

0.6185 

0.6655 

0.5356 

0.6185 

0.6705 

3 

0.2795 

0.6159 

0.5069 

0.5607 

0.6033 

0.4856 

0.5607 

*  0.6078 

4 

0.1668 

0.4663 

0.4365 

0.4829 

0.5196 

0.4182 

0.4829 

0.5234 

. 

5 

0.0853 

0.3831 

0.1751 

0.4149 

0.4464 

0.3593 

0.4149 

0.4497 

6 

0.0401 

0.4270 

0.3815 

0.4220 

0.4541 

0.3654 

0.4220 

0.4574 

7 

0.0051 

0.4507  ' 

0.3580 

0.3960 

0.4261 

0.3430 

0.3960 

0.4293 

8 

-0.0057 

0.2556 

0.2522 

0.2790 

0.3002 

0.2416 

0.2790 

0.3024 

9 

0.0084 

0.  1539 

0.1543 

0.1707 

0.1837 

0.1478 

0.1707 

0.1851 

to 

0.0375 

0.1195 

0.  1187 

0.1313 

0.1413 

0.1137 

0.1313 

0.1424 

11 

0.0604 

0.1325 

0.1165 

0.1289 

0.1387 

0.1116 

0.1289 

0.1397 

12 

0.0577 

0.1310 

0.1081 

0.1195 

0.1286 

0.1035 

0.1195 

0.1296 

13 

0.0475 

0.0836 

0.0951 

0.1052 

0.1132 

0.0911 

0.1052 

0.1140 

14 

0.0160 

0.1227 

0.0999 

0.1  105 

0.1189 

0.0957 

0.1105 

0.1198 

1? 

-0.0068 

0.1130 

0.0919 

0.1017 

0.1094 

0.0880 

0.1017 

0.1102 

16 

-0.0170 

0.0580 

0.0652 

0.0721 

0.0776 

0.0624 

0.0721 

0.0782 

IT 

-0.0260 

0.0594 

0.0531 

0.0538 

0.0588 
0.0595 

0.0633 
0.0640 

0.0509 
0.0515 

0.0588 
0.0595 

0.0637 
0.0645 

18 

-0.0325 

0.0583 

19 

-0.0214 

0.0620 

0.0515 

0.0569 

0.0613 

0.0493 

0.0569 

0.0617 

20 

-0.0027 

0.0455 

0.0431 

0.0477 

0.0514 

0.0413 

0.0477 

0.0517 

21 

0.0354 

0.0380 

0.0388 

0.0430 

0.0462 

0.0372 

0.0430 

0.0466 

22 

0.0533 

"  0.0504 

0.0418 

0.0463 

0.0498 

0.0401 

0.0463 

0.0502 

23 

0.0365 

0.0464 

0.0400 

0.0442 

0.0476 

0.0383 

0.0442 

0.0479 

24 

0.0242 

0.0337" 

"  0.0361 

0.0399 

0.0430 

0.0346 

0.0399 

0.0433 

25 

0.0232 

0.0458 

0.0366 

0.0405 

0.0436 

0.0351 

0.04C5 

0.0439 

26 

0.0122 

0.0365 

"  0.03)5 

0.0370 

0.0398 

0.0321 

0.0370 

0.0401 

27 

0.0137 

0.0292 

0.0268 

0.0296 

0,0318 

0.0256 

0.0296 

0.4)321 

28 

0.0032 

"  0.0235  " 

0.0250 

0.0276 

"0.0297 

0.0239 

0.0276 

0.0299 

29 

-0.0079 

0.0343 

0.0284 

0.0314 

0.0338 

0.0272 
0.0263 

0.0314 
'  0.0304 

0.0340 
0.0329 

30 

-0.0253" 

0^0335 

0.0275 

0.0304 

(5.032? 

31 

-0.0262 

0.0204 

0.0223 

0.0247 

0.0266 

0.0214 

0.0247 

0.0268 

32 

-0.0133 

0.0246 

0.0207 

0.0229 

"0.0246 

0.0196 

0.0229 

0.0248 

33 

-0.0067 

0.0221 

0.0201 

0.0222 
0.0209 

0.0239 
"0.0224 

0.0192 
0.0181 

0.0222 
0.0209 

0.0241 
0.0226  " 

34 

-0.0127 

0.0201 

0.0189 

35 

-0.0162 

0.0212 

0.0200 

0.0221 

0.0237 

0.0191 

0.0221 

0.0239 

-36" 

-0.0273 

6.0258 

0.0215 

0.0238 

0.0256 

0.0206 

0.0238 

0.0257 

37 

38" 

-0.0333 

0.0222 

0.0221 

0.0244 
0.0254  " 

0.0263 

"  o.02rj 

0.0211 
0.0220 

0.0244 
0.0254 

0.0265 
0.0275 

-0.0330 

0.0274 

0.022$ 

39 

-0.0391 

0.0245 

0.0209 

0.0232 

0,0249 

0.0201 

0.0232 

0.0251 

40 

-0.0413 

0.0163 

0.0184 

0.0204 

6.0219 

0.0176 

0.0204 

0.0221 

il8 


POWER  SPECTRUM  ESTIMATION 

FOR  OISCREIF  TIME  SERIES 

USING  TUKEY-MANNING  WEIGHTING  FACTORS 


TWO  SIZE  BEADS    40  ROWS  MIXED 


LAG 


40 


ttl 


SMOOTHED  POWER  SPECTRUM 
WITH  90  PCT.  CONFIDENCE  LIMITS 
12     3     4     5     6     7     8 


SMOOTHED  POWER  SPECTRUM 
KITH  95  PCT.  CONFIDENCE  LIMITS 
12     3     4     5     6     7 


0    1                                                                                                                          ••    ♦!         0    I                                                                                                                           .    •♦! 
11                                                                                                                       ••♦III                                                                                                                        .    •♦    I 

2  1                                                                                                                     .•♦121                                                                                                                   .       ♦          1 

3  1                                                                                                                .•♦131                                                                                                                .    •♦            1 

4  I                                                                                                             .•♦                 I4|                                                                                                  .         .    •♦                 I 

5  1                                                                                                     ••♦151                                                                                                      .    •♦                      I 

6  1                                                                                                        .*♦                      16    1                                                                                                      .       ♦                      I 
71                                                                                                     .•♦                         I         7    I                                                                                                      .♦! 

8  1                                                                                        .    •♦                                   18    1                                                                                         .    •♦                                   I 

9  1                                                                         .»♦                                                     19    1                                                                       .       ♦                                                     I 
10    I                                                               .    •♦                                                             1       10    1                                                                .    •♦                                                             I 

11  1                                                               ••♦                                                               I       11    1                                                               .    •♦                                                             1 

12  I                                                            .•♦                                                                  I       12    I                                                             .    •♦                                                               I 

13  1                                                       .    •♦                                                                    I       13    I                                                        .    •♦                                                                     1 

14  1                                                          ••♦                                                                    1       14    I                                                          .    •♦                                                                  I 

15  I                                                       .•♦                                                                       1       15    1                                                        .    ♦                                                                       I 

16  I                                           .    •♦                                                                                 I       16    I                                           .    •♦                                                                                 I 

17  1                                      .•♦                                                                                         f      17    1                                   .    •♦                                                                                         I 

18  I                                      .•♦                                                                                        t       18    I                                      .    ♦                                                                                         I 

19    1                                   .•    ♦                                                                                        I       19    I                                   .    •♦                                                                                         I 

?0    I                              ••♦                                                                                                I       20    I                              .    •♦                                                                                              I 

21  1                       .    •♦                                                                                             f      21    I                            •♦                                                                                              I 

22  I                          .    •♦                                                                                           I      22    1                          .    •♦                                                                                           I 

23  1                          .•♦                                                                                             1      23    1                          .    •*                                                                                           I 

24  I                      .    •♦                                                                                                     1       24    I                      .    •♦                                                                                                      I 

25  I                         .•♦                                                                                                     I      25    I                         •    ♦                                                                                                      I 

26  I                    .    •♦                                                                                                        I      26    I                    .    •♦                                                                                                        I 

27  r      ."•♦ r'zri .  •♦                                                                i 

28  1            .•♦                                                                                                                  I       28    I         .    •♦                                                                                                                   I 

29  I       '      .    •♦                                                                                                            1       29    1               .    •♦                                                                                                              I 

30  I              ••♦                                                                                                               1       30    I               .    •♦                                                                                                             1 

31  I       .•♦                                                                                                                       I      31    I       .    •♦                                                                                                                     1 

32  I    .•♦                                                                                                                          1       32    I    .    ■♦                                                                                                                        1 

»r.«" — -". "■ t"  3TTVV — - — — r 

34    !.•♦                                                                                                                            1       34    1.    •♦                                                                                                                          I 

35  I   .•♦        "~"    "~ ' " T "■■3¥"T"V*V " f 

36  1       ••♦                                                                                                                       1       36    I       ••♦                                                                                                                        I 

37  I       .•♦                                                                                                                       I       37    I       .    •♦                                                                                                                     I 

38  1       .    •♦                                                                                                                     I       38    1       .    •♦                                                                                                                     I 

jrT.w ■'" 1    39  I  .  •*    - • "T" 

40    !.•♦                                                                                                                     I      40    !.•♦                                                                                                                     I 

2   TIME 


1.3569  MIN 
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TWO  SIZE  BEAOS 


40  ROWS  MIXED 


LAG 


40 


SMOOTHED  SPECTRUM 

WITH  INDICATED  CONFIDENCE  LSPITS 

90  PCT. 

90  PCI. 

95  PCT. 

95  PCT. 

COVARIANCE 

RAW 

LOWER 

UPPER 

LOWER 

UPPER 

LAG 

SPECTRUM 

SPECTRUM 

LIMIT 

AVERAGE 

LIMIT 

LIMIT 

AVERAGE 

LIMIT 

0 

0.9748 

0.6457 

0.5549 

0.6138 

0.6605 

0.5316 

0.6138 

0.6654 

1 

0.7058 

0.5820 

0.5712 

0.6318 

0.6798 

0.5472 

0.6318 

0.6849 

2 

0.4705 

0.7177 

0.5446 

0.6025 

0.6483 

0.5217 

0.6025 

0.6531 

3 

0.2788 

0.3926 

0.4995 

0.5526 

0.5946 

0.4785 

0.5526 

0.5990 

4 

0.1308 

0.7075 

0.5338 

0.59C5 

0.6354 

0.5114 

0.5905 

0.6401 

5 

0.0277 

0.5545 

0.5319 

0.5884 

0.6331 

0.5095 

0.5884 

0.637B 

6 

-0.0383 

0.5369 

0.4693 

0.5191 

0.5585 

0.4495 

0.5191 

0.5627 

7 

-0.0681 

0.4480 

0.3783 

0.4184 

0.4502 

0.3624 

0.4184 

0.4536 

B 

-0.0766 

0.2408 

0.2665 

0.2948 

0.3172 

0.2553 

0.2948 

0.3195 

9 

-0.0500 

0.2495 

0.1885 

0.2085 

0.2243 

0.1805 

0.2085 

0.2260 

10 

-0.0153 

0.0941 

0.1307 

0.1446 

0.1556 

0.1252 

0.1446 

0.1567 

11 

0.0219 

0.1406 

0.1060 

0.1172 

0.1262 

0.1015 

0.1172 

0.1271 

12 

0.0435 

0.0937 

0.1035 

0.1  145 

0.1232 

0.0991 

0.1145 

0.1241 

13 

0.0376 

0.1299 

0.0986 

0.1090 

0.1173 

0.0944 

0.1090 

0.1182 

1* 

0.0379 

0.0826 

0.0871 

0.0963 

0.1036 

0.0834 

0.0963 

0.1044 

15 

0.0214 

0.0902 

0.0745 

0.0824 

0.0887 

0.0714 

0.0824 

0.0894 

16 

0.0149 

0.0668 

0.0619 

0.0685 

0.0737 

0.0593 

0.0685 

0.0742 

17 

0.0164 

0.0501 

0.0501 

0.0555 

0.0597 

0.0480 

0.0555 

0.0601 

IB 

0.0054 

" 0.0548 

"  0.0483 

0.0534 

0.0574 

0.0462 

0.0534 

0.0579 

19 

0.0020 

0.0537 

0.0447 

0.0495 

0.0532 

0.0428 

0.0495 

0.0536 

20 

-0.0065 

0.0356 

0.0395 

0.0437 

0.0470 

0.0378 

0.0437 

0.0474 

21 

-0.01C5 

0.0500 

0.0384 

0.0424 

0.045f 

0.0368 

0.0424 

0.0460 

22 

-0.0158 

0.0343 

0.0369 

0.0408 

0.0439 

0.0354 

0.0408 

0.0443 

23 

-0.0268 

0.0448 
0.0333 

0.0355 
0.0322 

0.0393 
0.0356 

0.0423 
0.0383 

0.0340 
0.0308 

0.0393 
0.0356 

0.0426 

\      2* 

-0.0227 

0.0386 

i   25 

-0.0250 

0.0309 

0.0280 

0.0310 

0.0333 

0.0268 

0.0310 

0.0336 

""26 

-0.0252 

0.0289 

0.0281 

0.0311 

0.0334 

0.0269 

0.0311 

0.0337 

1  27 

-0.0305 

0.0358 

0.0283 

0.0313 

0.0337 

0.0271 

0.0313 

0.0340 

T28 

-0.0466 

0.0250 

0.0275 

0.0304 

0.0327 

0.0263 

0.0304 

0.0329 

L   2" 

-0.0500 
-0.0460 

0.0358 

0.0262 

0.0290 
0.0242 

0.0312 
0.0261 

0.0251 
0.0210 

0.0290 
0.0242 

0.0315 

30 

0.0196 

0.0219 

0.0262 

31 

-0.0388 

0.0219 

0.0180 

0.0199 

0.0215 

0.0173 

0.0199 

0.0216 

32 

-0.0385 

0.0164  "' 

0.0204 

0.0226 

0.0243 

0.0196 

0.0226 

0.0245 

33 

-0.0142 

0.0358 

0.0261 

0.0289 

0.0311 

0.0250 

0.0289 

0.0313 

34 

0.0189 

0.0276 

0.0261 

0.0289 

0.0311 

0.0250 

0.0289 

0.0313 

35 

0.0533 

0.0246 

0.0223 

0.0247 

0.0266 

0.0214 

0.0247 

0.0268 

36 

0.0777 

0.0220" 

0.0208 

0.0230 

0.0247  ~ 

"   0.0199 

0.0230 

0.0249 

37 

0.0805 

.  0.0232 

0.0194 

0.0215 

0.0231 

0.0186 

0.0215 

0.0233 

38 

0.0726 

"0.0175 

0.0169 

0.0210 

0.0225 

0.0181 

0.0210 

0.0227 

39 

0.0557 

0.0256 

0.0194 

0.0215 

0.0231 

0.0186 

0.0215 

0.0233 

40 

0.0262 

0.0171 

0.0193 

0.0214 

0.0230 

0.0185 

0.0214 

0.0232 

• 
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TWO  SIZE  BEADS    *0  ROWS  MIXFP 


LAG  -  '.0 


C*) 


SMCOTHED  POWFR  SPECTRUM 
WITH  9C  PCT.  CCNriDENCE  LIMIIS 
2     J     4     5     ft     7     R 


0 

1 

I 
2 
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3 

4 

5 

6 

T 

I 

8 

9 

10 

11 

12 

13 

14 

15 

16 

1/ 

18 

19 

20 

21 

22 

23 

24 

[        . 

25 

I       .  •  ♦ 

26 

1       •  •♦ 

27 

1        •  M 

28 

1      .  M 

29 

1      •  •♦ 

™ 

[   .  •♦ 

31 

!••♦ 

32 

1  .  •♦ 

33 

1      .•♦ 

34 

1       »•♦ 

35 

1   .  •♦ 

36 

I  •  •  ♦ 

37 

1  ••♦ 

3B 

1.  •♦ 

39 

I  .•♦ 

40 

!  .  *4 

•  •♦ 
.44 

.  •♦ 
•  4 


SMCOTHFO    PCWER     SPECTRUM 
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12  3  4  5  6  7  8 
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I    3 
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1 
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1 

1   14 

1 
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1 

1   16 

1 

. 

1   17 

1 

.    4 

1   18 

1 

.  M 

1   19 

1 

.  »4 
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1 

. 

•  4 

I   21 

1 

. 

•  4 

1   22 

1 

. 

•  4 

1   23 

1 

. 

■  4 

1   24 

1 

•  •  ♦ 

1   25 

I 

4 

I   26 

1       , 

♦ 

I   27 

1 

•  4 

I   2R 

1 

•  4 

1   29 

1 

•  4 

1   30 

1   .  •♦ 

1   31 

!•■♦ 

1   32 

I  .  •♦ 

I   33 

1 

4 

1   34 

1 

•  4 

I   35 

1   .  •♦ 

1   36 

I  .  •♦ 

1   37 

1  .  ♦ 

1   38 

1.  •♦ 

1   39 

I.  •♦ 

1   40 

1.  •♦ 

.       »4| 
.      «4 


.       •♦ 

•  4 

•  4 
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IHREE    Sl/F     PFACS 

to 


*o  rows  Mtxrn 


I  40 


SMUnfMEl'  SPECTRUM 

WITH  IND 

ICATEO  CUNFIDENCF  L  1  f  1  TS 

.  . 

90  PCI . 

90  PCI. 

95  PCI. 

95  PCI 

COVARIANCE 

RAM 

LC-I  R 

UPPER 

LUwbR 

UPPER 

LAG 

SPFG  TRUM 

SPECTRUM 

LIMII 

Avr RAGE 

I  |MI  1 

LIM  1 

AVERAGE 

LIMII 

0 

0.9035 

0.4561 

0.43  M 

0.48  15 

0.5203 

0.4187 

0.4835 

0.5242 

1 

0.6015 

0.5110 

0.4  1*7 

0.48C3 

0.5168 

0.4159 

0.4803 

0.5206 

? 

0.  3676 

0.*430 

0.4  )67 

I) .  '.  8  3  1 

0.5198 

0.4183 

C.4B31 

0.5237 

3 

0.2012 

0.5353 

0.4675 

0.S17? 

0.5565 

0.44  79 

0.517? 

0.5606 

* 

0.0977 

0.5550 

C.5005 

1.S516 

0.5957 

0.479* 

0.5536 

0.6001 

5 

0.0?*t 

0.5691 

0.4802 

U  .  5  1 1  7 

0.5715 

0.4600 

0.5112 

0.5758 

6 

-0.0311 

0.431  3 

0.  3913 

0.4  194 

0.4728 

o.3eot 

0.4394 

0.4764 

7 

-0.0689 

0. 1760 

0.2918 

0.1277 

0. 3473 

0.7795 

0.3227 

0. 3499 

8 

-O.0B16 

0.?076 

0.71 11 

0.7357 

0.7516 

0.20*1 

0.2357 

0.2555 

9 

-0.07*? 

0.7015 

n.l 719 

0. 1902 

0.2046 

0.16*7 

0.  1902 

0.2062 

10 

-0.06*0 

0. 1502 

0  .  1  '.  1  7 

0. 1 56  8 

0.168  7 

0.1358 

0.1568 

0. 1699 

It 

-0.0290 

0.1251 

0.1219 

0.1  148 

0.1451 

0.1168 

0.1348 

0. 1461 

I? 

-0.00*5 

0. 1388 

0.1261 

0.1 315 

0. 1502 

0.1208 

0.1395 

0.1513 

13 

0.0155 

0.1554 

o.i ?in 

0.1  338 

0. 14*0 

0.1  159 

0.1338 

0. 1451 

1* 

0.0291 

0.0857 

0 .0929 

0.102R 

0.1106 

0.089C 

0.1028 

0.1114 

IS 

0.0215 

0.0844 

O.n  ?*4 

0.0M45 

0.0909 

0.0732 

0.0845 

0.0916 

16 

0.0150 

0.0815 

0.07*7 

0.0876 

0.0889 

0.0715 

0.0826 

0.0896 

17 

0.021? 

0.0791    ■ 

0.0776 

0.1H03 

0.0864 

0.0696 

0.08C3 

0.0871 

IS 

0.0311 

0.0797 

. 0.0652 

0.0712 

0.0788 

0.061* 

0.0732 

0.079* 

19 

0.0222 

0.0546 

0.0545 

0.0603 

0.0649 

0.0523 

O.0601 

0.0654 

20 

0.0099 

0.05?6 

0.04  17 

0.0528 

0.0568 

0.0*57 

0.0528 

0.0572 

?1 

-0.0087 

0.051 3 

0.0465 

0.0515 

0.0554 

0.0**6 

0.0515 

0.0558 

22 

-0.0122 

0.0 508 

0.O440 

0.0487 

0.0524 

0.0*22 

0.0487 

0.0528 

23 

-0.00*4 

0.0420 

0.0  386 

0.0477 

0.0459 

0.037C 

0.0427 

0.0463 

2* 

0.0009 

0.0361 

0.0  ISO 

0.0387 

0.0416 

0.0335 

0.0387 

0.0419 

25 

0.0088 

0.0406 

0.0376 

0.0361 

0.0388 

0.0312 

0.0361 

0.0391 

?6 

0.01C3 

0.0271 

0.0315 

0.0348 

0.0375 

0.0302 

0.0348 

0.0378 

27 

0.0014 

0.0446 

0.0353 

0.0190 

0.0420 

0.0338 

0.0390 

0.0423 

28 

-0.0139 

0.0397 

0.0354 

0.0391 

0.0421 

0.0339 

0.0391 

0.0424 

29 

-0.0078 

0.0325 

0.0309 

0.0342 

0.0368 

0.0296 

0.034? 

0.0371 

30 

-0.0062 

0.0320 

0.0280 

0.0310 

0.0333 

0.0268 

0.0310 

0.0336 

31 

-0.0009 

0.0273 

0.0?56 

0.0283 

0.0304 

0.0245 

0.0283 

0.0306 

32 

0.0053 

0.0264 

0.0246 

0.027? 

0.0293 

0.0236 

0.0272 

0.0295 

33 

-0.0058 

0.0287 

0.0243 

0.O269 

0.0289 

0.0233 

0.0269 

0.0291 

3* 

-0.0093 

0.0237 

0.0247 

0.0273 

0.0294 

0.0236 

0.0273 

0.0296 

35 

-0.0027 

0.0331 

0.0273 

0.03C2 

0.0325 

0.026? 

0.0302 

0.0327 

36 

-0.0024 

0.0310 

0.0272 

0.0301 

0.0324 

0.0261 

0.0301 

0.0326 

37 

-0.0123 

0.0253 

0.0240 

0.0266 

0.0286 

0.0230 

0.0266 

0.0?88 

38 

-0.0208 

0.0247 

0.0228 

0.0252 

0.0271 

0.0218 

0.025? 

0.0?73 

39 

-0.0281 

0.0261 

0.022B 

0.0253 

0.0272 

0.0219 

0.0253 

0.0274 

40 

-0.0367 

0.0241 

0.0227 

0.0251 

0.0270 

0.021T 

0.0251 

0.0272 
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THREE  SI2E  BEADS 

to 


40  ROUS  MIXFD 


LAC 


SMOOTHED  PCWER  SPECTRUM 
WITH  90  PCT.  CONFIDENCE  LIMITS 
I     2     3     *     5     6     7 


SMCOTHEO  POWER  SPECTRUM 
WITH  95  PCI.  CONFIDENCE  LIMITS 
12     3     4     5     6     7     9 


0  I 

1  1 

1  1 

2  1 

3  1 

*  1 

5  I 

6  I 

T  1 

8  1 

9  1 
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13  1 

14  1 

15  I 

16  1 

17  I 
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20  I 

21  1 

22  I 
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24  I 

25  I 

26  1 

. 

27  I 

28  I 

29  I 

. 

30  I 

.* 

31  1 

.  •♦ 

32  1 

•  •♦ 

33  I 

•  •♦ 

34  1 

•  ■  ♦ 

35  I 

»  • 

36  I 

.  • 

37  t. 

•  ♦ 

38  1. 

•  ♦ 

39  I. 

•  ♦ 

40  I. 

•  ♦ 

•  •  • 


•  •  ♦ 

.  •  ♦ 


1    0  1 

1     2  1 

1    3  1 

1     4  1 

1    5  I 

1    6  1 

1    7  1 

I    P  1 

1    9  1 
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1    III 

1   12  1 

1   1)1 

1    15  1 

1    16  1 

1    17  1 

1    IH  1 

1    19  1 

1    20  1 

I    21   1 

1,   22  I 

I   23  1 

1   24  1 

I   25  1 

1   26  1 

1   27  1 

I   28  1 

I   29  1 

1   30  1 
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1   32  1 

1   33  1 

1   34  I 

1   35  1 

1   36  1 

1   37  I. 

1   38  I. 

1   39  I. 

1   40  t. 
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SMrniHFO    SPECTRI/M    WITH     INDICATED    CONFIDENCE    UMTS 


90  PCT  . 

90  PCT. 

95  PCT. 

95  PCT. 

COVAR (ANCE 

■    RAW 

u:wfr 

UPPER 

LOWER 

UPPER 

LAG 

SPFC  TRUM 

SPFC TRUM 

LIMIT 

avfrage 

LIMIT 

LIMIT 

AVERAGE 

LIMIT 

0 

0.9281 

0.4511 

0.4741 

0.4694 

0.5050  , 

0.4C65 

0.4694 

0.5088 

1 

0.6461 

0.4876 

0.4451 

0.4  9  74 

0.5298 

0.4264 

0.4924 

0.5338 

2 

0.4014 

0.5434 

0.4  865 

0.5  IP? 

0.5791 

0.466C 

0.53e? 

0.5834 

3 

0.220* 

0.5783 

0.504? 

(1.5  5  78 

0.6002 

0.4830 

0.5578 

0.6046 

S 

0.0973 

0.5311 

0.5717 

U  .  5  7  7  1 

0.6209 

0.4997 

f) .  5  7  1 1 

0.6255 

5 

0.01  JO 

0.6677 

0.531? 

0.5876 

0.6373 

O.5089 

0.5876 

0.6370 

6 

-0.0351 

0.4  8  39 

0.4  150 

0.481? 

0.51 78 

0.4168 

0.4812 

■  0.5217 

7 

-0.0711 

0.2895 

0.787'. 

J.  1179 

0. 3471 

0.2753 

0.31 79 

0.3446 

8 

-0.0863 

0.2087 

0.70?'. 

0.7739 

0.2409 

0.1939 

0.2239 

0.2427 

9 

-0.0R18 

0. 1884 

0.1824 

0.7018 

0.71 71 

0.1 74e 

0.2018 

0.2188 

10 

-0.0773 

0.2216 

0.1781 

0. 1970 

0.21  19 

0.17C6 

0.1970 

0.2135 

11 

-0.0682 

0. 156? 

0.155? 

0.1717 

0.1848 

0.1487 

0.1717 

0.1861 

1Z 

-0.03*1 

0.1578 

0.128? 

0.1418 

0. 1526 

0.1228 

0. 141R 

0. 1537 

13 

0.0101 

0.1055 

0.1005 

0.1111 

0.1 196 

0.0962 

0.1111 

0.1205. 

1* 

0.0305 

0.0808 

0.0780 

0.0863 

0.0929 

0.074P 

0.0863 

0.0936 

15 

0.0496 

0.0781 

0.0719 

O.ORIR 

0.08R0 

0.0708 

0.0818 

0.08H7 

16 

0.0374 

0.0898 

0.0754 

0.OR34 

0.0897 

0.0722 

0.0834 

0.0904 

17 

0.0146 

0.0756 

0.0701 

0.0776 

0.0835 

0.0672 

0.0776 

0.0841 

IB 

O.O0C6 

0.069J 

0.0586 

0.064R 

0.0697 

0.0561 

0.0648 

0.0703 

19 

-0.0117 

0.0450 

0.048'. 

0.0536 

0.0577 

0.0464 

0.0536 

0.0581 

70 

-0.0719 

0.0550 

0.0445 

0.0492 

0.0529 

0.0426 

0.0492 

0.0533 

21 

-0.04C9 

0.0418 

0.04C  1 

0.0446 

0.0480 

0.0386 

0.0446 

0. 04  8  3 

?2 

-0.0361 

0.0398 

0.0336 

0.0372 

0.0400 

0.0322 

0.037? 

0.040  3 

23 

-0.0300 

0.0273 

0.0317 

0.0350 

0.0377 

0.0303 

0.0350 

0.0380 

2*. 

-0.0171 

0.0458 

0.0359 

0.0397 

0.0427 

0.0344 

0.0397 

0.0430 

25 

-0.0093 

0.0399 

0.0375 

0.0415 

0.0446 

0.0359 

0.0415 

0.0449 

26 

-0.0001 

0.0402 

0.0347 

0.03P4 

0.0413 

0.0333 

0.03e4 

0.0416 

27 

0.0065 

0.0333 

0.0300 

0.0332 

0.0357 

0.0287 

0.0332 

0.0360 

2ft 

0.0055 

0.0258 

0.0747 

O.0?73 

0.0794 

0.0236 

0.0273 

0.0296 

29 

0.0096 

0.0241 

0.0230 

0.0255 

0.0774 

0.0221 

0.0255 

0.0276 

30 

0.0225 

0.0278 

0.0231 

0.0256 

0.0275 

0.0222 

0.0256 

0.0277 

31 

0.0228 

0.0226 

0.0738 

0.0263 

0.0283 

0.0228 

0.0263 

0.0285 

3? 

0.0294 

0.0324 

0.0249 

0.0275 

0.0296 

0.0238 

0.0275 

0.0298 

33 

0.0209 

0.0227 

0.0214 

0.0237 

0.0255 

0.0205 

0.0237 

0.0257 

3* 

0.0073 

0.0172 

0.0177 

0.0196 

0.0211 

0.0170 

0.0196 

0.0213 

35 

-0.0038 

0.0214 

0.0195 

0.0216 

0.0232 

0.0187 

0.0216 

0.0234 

96 

-0.0111 

0.0263 

0.0231. 

0.0255 

0.0275 

0.0221 

0.0255 

0.0277 

37 

-0.0184 

0.0280 

0.0238 

0.0264 

0.0284 

0.0228 

0.0264 

0.0286 

36 

-0.0232 

0.0232 

0.0214 

0.0237 

0.0255 

0.0205 

0.0237 

0.0257 

39 

-0.0141 

0.0204 

0.0209 

0.0232 

0.0249 

0.0201 

0.0232 

0.0251 

♦0 

O.O002 

0.0287 

0.0222 

0.0246 

0.0264 

0.0213 

0.0246 

0.0266 
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90  PCI 

C0VAR1ANCF 

RAW 

LOWER 

LAG 

SPECTRUM 

SPECTRUM 

LIMIT 

0 

0.8960 

0.4985 

0.*I2P 

0.4939 

0.41*8 

0. 3T48 

0.2325 

0.4190 

0.3541 

0.0761 

0.3141 

0.  1003 

0.0184 

0.3215 

0.3090 

'COO*") 

0.4102 

0.3251 

0.011* 

0.2966 

0.2851 

0.0276 

0.25R* 

0.2473 

0.028T 

0.2810 

0.2475 

0.0223 

0  .  ?  74  'J 

0.2460 

0.0125 

0.2578 

0.2307 

0.0036 
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SMGnTMrn  SPECTRUM 

WITH  INDICATED  CONFIDENCE  LIMITS 

90  PCT. 

90  PCT. 

95  PCT. 

95  PCT. 

COVARIANCE 

RAW 

LOWER 

UPPER 

LOWER 

UPPER 

LAG 

SPECTRUM 

SPECTRUM 

LIMIT 

AVERAGE 

LIMIT 

LIMIT 

AVERAGE 

LIMIT 

0 

0.9378 

0.3632 

0.3367 

0.3775 

0.4008 

0.3225 

0.3725 

0.4037 

0.5003 

0.381 7 

0.3256 

0.3602 

0.3876 

0.3119 

0.3602 

0.3904 

0.2076 

0.314  1 

0.3395 

0.3755 

0.4041 

0.3252 

0.3755 

0.4071 

0.0160 

0.4473 

0.36S7 

0.4045 

0.4352 

0.3503 

0.4045 

0.4385 

-0.0542 

0.3193 

0.3115 

0.3446 

0.3708 

0.2984 

0.3446 
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0.7570 
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R 
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0.3108 

•   0.2502 

0.7889 

0.3132 

9 

0.0265 

0.3115 

0.7886 

0.3193 

0.3435 

0.2765 

0.3193 

0.3461 

10 

-0.0284 

0. 3686 

0.3366 

0.3773 

0.4006 

0.3224 

0.  1723 

0.4036 

11 

-0.0670 

0.4406 

0.3343 

0.3698 

0.3979 

0.3202 

0.  1698 

0.4009 

12 

-0.063* 

0.7793 

0.7334 

0.75P? 

0.2778 

0.2236 

0.7582 
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I) 

-0.0356 

0.1336 
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0.1417 

0.1636 

0.17  74 
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0.1746 

0.1 3  78 

0.1483 

0.1 194 

0. 1 378 
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15 
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0. 104S 

0.1 156 

0.1244 

0.1C01 

0. 1156 

0.1253 

16 

0.0163 

0. 1013 

0.0976 

0. 1075 

0.1103 

0.0887 

0.1025 

O.lllt 

17 

-0.0261 

0.1057 

0.0835 

0.0973 

0.0994 

0.0800 

0.0923 

0.1001 

IB 

-0.0426 

0.0566 

0.0692 

0.0765 

0.0824 

0.0663 

0.0765 

0.0830 

19 

-0.0340 

0.0872 

0.0739 

0.0817 

0.0879 

0.07C8 

0.0817 

0.0886 

20 

-0.0078 

0.0959 

0.0744 

0.0873 

0.0886 

0.0713 

0.0823 

0.0892 

21 

0.0252 

0.0503 

0.0660 

0.0731 

0.0786 

0.0633 

0.0731 

0.0792 

72 

0.0457 

0.0958 

0.0636 

0.0703 

0.0757 

0.0609 

0.0703 

0.0762 

21 

0.0485 

0.0394 

0.0544 

0.0602 

0.0648 

0.0521 

0.0602 

0.0653 

24 

0.0337 

0.0662 

0.0529 

0.0585 

0.0629 

0.0506 

0.0585 

0.0634 

25 

0.0168 

0.0620 

0.0552 

0.0610 

0.0657 

0.0529 

0.0610 

0.0662 

26 

0.0096 

0.0538 

0.04B1 

0.0532 

0.0573 

0.0461 

0.0532 

0.0577 

27 

0.0069 

0.0432 

0.0431 

0.0477 

0.0513 

0.0413 

0.0477 

0.0517 

28 

0.0445 

0.0504 

0.0426 

0.0471 

0.0507 

0.0408 

0.0471 

0.0511 

29 

0.0595 

0.0446 

0.0429 

0.0474 

0.0510 

0.0411 

0.0474 

0.0514 

10 

0.0164 

0.0502 

0.045R 

0.0506 

0.0545 

0.0438 
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0.0549 
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0.0497 
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0.0487 
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3* 
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0.0470 

0.0506 

0.0407 

0.0470 

0.0510 
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-0.0295 

0.0429 

0.0411 

0.0455 

0.0489 

0.0394 

0.0455 

0.0493 

36 

-0.0031 

0.0462 

0.0386 

0.0427 

0.0459 

0.0369 

0.0427 

0.0462 

37 

-0.0134 

0.0354 

0.0367 

0.0405 

0.0436 

0.0351 

0.0405 

0.0440 
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-0.0295 

0.0452 

0.0369 

0.0408 

0.0439 

0.0354 
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0.0443 

39 

-0.0343 
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0.0398 
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0.0474 

0.0382 
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